Minor Project: Oocyte Alternative Splicing Analysis

1 Introduction to this Analysis

1.1 Part 1: Upstream Analysis

The upstream analysis is performed on the CRG cluster using the Bash command line. This phase involves:

1. Data Acquisition:
   • Downloading *.fastq.gz files from 3 different databases.
2. Alignment:
   • Aligning the sequences to the latest mouse genome release (mm10) using Bowtie2 through the Vast-Tools align function.
3. Event Detection:
   • Importing inclusion tables generated by Vast-Tools into RStudio. The tables cover the following splicing events:
     • Exons (EX)
     • Introns (IN)
     • Alternative 5’ Splice Sites (Alt5)
     • Alternative 3’ Splice Sites (Alt3)
     • Microexons (MIC)

---

1.2 Part 2: Statistical Analysis

1.2.1 Objective

To identify splicing events that are differentially spliced in Oocytes treated with PladB with two differentes doses, or modified to be FMN2 DKO or SPIRE1/2 DKO.

1.2.2 Steps

1. Statistical Testing:
   • Use the betAS package in R to perform:
     • False Discovery Rate (FDR)
     • F-statistic
     • Probability of Differential Splicing (Pdiff)
2. Filtering Criteria:
   • Retain events with:
     • FDR ≤ 0.05
     • Pdiff ≥ 0.95
3. Batch Effect Mitigation:
   • A final list is generated by intersecting the splicing events identified in all five databases. This ensures the events are associated exclusively with differentiation. All exon related events (C1,C2,C3,S,including microexons) will be collapsed into a single list, EX.

---

1.3 Importing Inclusion Data

pladb_data <- getDataset(pathTables = paste0(getwd(),"/inclusion_tables/pladb_INCLUSION_LEVELS_FULL-mm10.tab"), tool = "vast-tools")
pladb_events <- getEvents(pladb_data, tool = "vast-tools") # Extract alternative splicing events
pladb_events <- alternativeEvents(pladb_events, minPsi = 1, maxPsi = 99) # Filter events based on PSI thresholds
pladb_exons <- filterEvents(pladb_events, types = c("C1", "C2", "C3", "S", "MIC"), N = 3)
pladb_alt <- filterEvents(pladb_events, types = c("Alt5", "Alt3"), N = 3)
pladb_introns <- filterEvents(pladb_events, types = c("IR"), N = 3)


fmndko_data <- getDataset(pathTables = paste0(getwd(),"/inclusion_tables/fmndko_INCLUSION_LEVELS_FULL-mm10.tab"), tool = "vast-tools")
fmndko_events <- getEvents(fmndko_data, tool = "vast-tools") # Extract alternative splicing events
fmndko_events <- alternativeEvents(fmndko_events, minPsi = 1, maxPsi = 99) # Filter events based on PSI thresholds
fmndko_exons <- filterEvents(fmndko_events, types = c("C1", "C2", "C3", "S", "MIC"), N = 3)
fmndko_introns <- filterEvents(fmndko_events, types = c("IR"), N = 5)
fmndko_alt <- filterEvents(fmndko_events, types = c("Alt5", "Alt3"), N = 3)



spire_data <- getDataset(pathTables = paste0(getwd(),"/inclusion_tables/spiredko_INCLUSION_LEVELS_FULL-mm10.tab"), tool = "vast-tools")
spire_events <- getEvents(spire_data, tool = "vast-tools") # Extract alternative splicing events
spire_events <- alternativeEvents(spire_events, minPsi = 1, maxPsi = 99) # Filter events based on PSI thresholds
spire_exons <- filterEvents(spire_events, types = c("C1", "C2", "C3", "S", "MIC"), N = 3)
spire_introns <- filterEvents(spire_events, types = c("IR"), N = 5)
spire_alt <- filterEvents(spire_events, types = c("Alt5", "Alt3"), N = 3)

1.4 Metadata

1.4.1 PladB

# Load metadata file containing sample information
metadata_pladb <- read.csv(paste0(getwd(),"/metadata/metadata_pladb.csv"), sep = "\t")
metadata_pladb<-metadata_pladb[seq(1, nrow(metadata_pladb), by = 2),]
DT::datatable(metadata_pladb, options = list(pageLength = nrow(metadata_pladb), scrollX = TRUE))

1.4.2 FMN2DKO

# Load metadata file containing sample information
metadata_fmn2dko<-data.frame(samples=c("SRR6026682_SRR6026683_SRR6026684_merged_trimmed","SRR6026685_SRR6026686_SRR6026687_merged_trimmed","SRR6026688_SRR6026689_SRR6026690_merged_trimmed","SRR6026691_SRR6026692_SRR6026693_merged_trimmed"),condition=c("control","control","fmndko","fmndko"))
DT::datatable(metadata_fmn2dko, options = list(pageLength = nrow(metadata_fmn2dko), scrollX = TRUE))

1.4.3 Spire

# Load metadata file containing sample information
metadata_spire<-data.frame(samples=c("Oocytes_FG_Spire12_Cont_a","Oocytes_FG_Spire12_Cont_b","Oocytes_FG_Spire12_Cont_c","Oocytes_FG_Spire12_DKO_a","Oocytes_FG_Spire12_DKO_b","Oocytes_FG_Spire12_DKO_c"), condition=rep(c("control", "spiredko"), each=3))
DT::datatable(metadata_spire, options = list(pageLength = nrow(metadata_spire), scrollX = TRUE))

---

2 Splicing & Sample Quality Checking

Splicing inclusion plots look like a U when looking at the exons, as normally most of them are spliced out completely (dPSI=0) or constitutively spliced in (dPSI=100), but some of them are in between and those are the most interesting for the analysis. On the other hand, introns are ussually not included (that is their definiiton), but in a few cases they are included. Thus, the intron plot looks like the exponential distribution.

2.1 PSI Distribution Plots

2.1.1 PladB

# Create the first plot for exon inclusion levels
bigPicturePlotExons <- bigPicturePlot(table = pladb_exons$PSI)

# Customize the first plot for professional reporting
plot1 <- bigPicturePlotExons +
  ggtitle("Exon Inclusion Levels") +
  xlab("Samples") +
  ylab("PSI Values") +
  theme_minimal() +
  theme(
    plot.title = element_text(hjust = 0.5, face = "bold", size = 14),
    axis.text = element_text(size = 10),
    axis.title = element_text(size = 12),
    panel.background = element_rect(fill = "lightblue"),  # Set background color
    panel.grid.major = element_line(color = "gray90", size = 0.5),
    panel.grid.minor = element_line(color = "gray95", size = 0.25)
  )

# Create a second plot for a hypothetical dataframe
# Replace 'hypothetical_data$PSI' with your actual dataframe and column
bigPictureIntrons <- bigPicturePlot(table = pladb_introns$PSI)

# Customize the second plot
plot2 <- bigPictureIntrons +
  ggtitle("Introns Inclusion Levels") +
  xlab("Samples") +
  ylab("PSI Values") +
  theme_minimal() +
  theme(
    plot.title = element_text(hjust = 0.5, face = "bold", size = 14),
    axis.text = element_text(size = 10),
    axis.title = element_text(size = 12),
    panel.background = element_rect(fill = "lightblue"),  # Set background color
    panel.grid.major = element_line(color = "gray90", size = 0.5),
    panel.grid.minor = element_line(color = "gray95", size = 0.25)
  )

# Combine both plots into one figure using patchwork
plot1 / plot2

2.1.2 Spire

# Create the first plot for exon inclusion levels
bigPicturePlotExons <- bigPicturePlot(table = spire_exons$PSI)

# Customize the first plot for professional reporting
plot1 <- bigPicturePlotExons +
  ggtitle("Exon Inclusion Levels") +
  xlab("Samples") +
  ylab("PSI Values") +
  theme_minimal() +
  theme(
    plot.title = element_text(hjust = 0.5, face = "bold", size = 14),
    axis.text = element_text(size = 10),
    axis.title = element_text(size = 12),
    panel.background = element_rect(fill = "lightblue"),  # Set background color
    panel.grid.major = element_line(color = "gray90", size = 0.5),
    panel.grid.minor = element_line(color = "gray95", size = 0.25)
  )

# Create a second plot for a hypothetical dataframe
# Replace 'hypothetical_data$PSI' with your actual dataframe and column
bigPictureIntrons <- bigPicturePlot(table = spire_introns$PSI)

# Customize the second plot
plot2 <- bigPictureIntrons +
  ggtitle("Introns Inclusion Levels") +
  xlab("Samples") +
  ylab("PSI Values") +
  theme_minimal() +
  theme(
    plot.title = element_text(hjust = 0.5, face = "bold", size = 14),
    axis.text = element_text(size = 10),
    axis.title = element_text(size = 12),
    panel.background = element_rect(fill = "lightblue"),  # Set background color
    panel.grid.major = element_line(color = "gray90", size = 0.5),
    panel.grid.minor = element_line(color = "gray95", size = 0.25)
  )

# Combine both plots into one figure using patchwork
plot1 / plot2

2.1.3 FMN2DKO

# Create the first plot for exon inclusion levels
bigPicturePlotExons <- bigPicturePlot(table = fmndko_exons$PSI)

# Customize the first plot for professional reporting
plot1 <- bigPicturePlotExons +
  ggtitle("Exon Inclusion Levels") +
  xlab("Samples") +
  ylab("PSI Values") +
  theme_minimal() +
  theme(
    plot.title = element_text(hjust = 0.5, face = "bold", size = 14),
    axis.text = element_text(size = 10),
    axis.title = element_text(size = 12),
    panel.background = element_rect(fill = "lightblue"),  # Set background color
    panel.grid.major = element_line(color = "gray90", size = 0.5),
    panel.grid.minor = element_line(color = "gray95", size = 0.25)
  )

# Create a second plot for a hypothetical dataframe
# Replace 'hypothetical_data$PSI' with your actual dataframe and column
bigPictureIntrons <- bigPicturePlot(table = fmndko_introns$PSI)

# Customize the second plot
plot2 <- bigPictureIntrons +
  ggtitle("Introns Inclusion Levels") +
  xlab("Samples") +
  ylab("PSI Values") +
  theme_minimal() +
  theme(
    plot.title = element_text(hjust = 0.5, face = "bold", size = 14),
    axis.text = element_text(size = 10),
    axis.title = element_text(size = 12),
    panel.background = element_rect(fill = "lightblue"),  # Set background color
    panel.grid.major = element_line(color = "gray90", size = 0.5),
    panel.grid.minor = element_line(color = "gray95", size = 0.25)
  )

# Combine both plots into one figure using patchwork
plot1 / plot2

2.2 Total Splicing Events

As it can be observed below, the dominic dataset has the highest number of splicing events, what can be explained by the higher number of reads (90M in myotubes and 140M in myoblasts).

# Adjusted version of the script
splicing_events <- tibble(
  event_type = names(spire_events$EventsPerType),
  pladb = pladb_events$EventsPerType,
  spire = spire_events$EventsPerType,
  fmndko = fmndko_events$EventsPerType
)

splicing_events_long <- splicing_events %>%
  pivot_longer(cols = c(pladb, spire, fmndko),
               names_to = "study", values_to = "count") %>%
  mutate(event_type = ifelse(event_type %in% c("C1", "C2", "C3", "ANN", "S", "MIC"), "EX", event_type)) %>%
  group_by(study, event_type) %>%
  summarize(count = sum(count, na.rm = TRUE), .groups = "drop")

# Proportion bar plot
plot_proportion <- ggplot(splicing_events_long, aes(fill = event_type, y = count, x = study)) +
  geom_bar(position = "fill", stat = "identity") +
  labs(
    title = "<b style='font-size:18px;'>Proportion of Splicing Events by Study</b>",
    x = NULL,
    y = "Proportion of Events",
    fill = "Event Type"
  ) +
  theme_minimal(base_family = font) +
  theme(
    axis.text.x = element_blank(),
    axis.ticks.x = element_blank(),
    plot.title = element_markdown(),
    panel.background = element_rect(fill = "gray98"),
    panel.grid = element_line(color = "gray90", size = 0.5)
  ) +
  scale_fill_paletteer_d("MetBrewer::Hokusai3")

splicing_events_long_sum<-splicing_events_long %>%
  dplyr::group_by(study) %>%
  summarize(count = sum(count))


# Absolute counts plot
plot_absolute <- ggplot(splicing_events_long_sum, aes(y = count, x = study)) +
  geom_point(color = "#C70039", size = 3) +
  geom_text(aes(label = count), vjust = -1, size = 4, check_overlap = TRUE) +
  labs(
    x = "Study",
    y = "Number of Events",
    caption = paste0("Created by AG on ", Sys.Date()), size=3
  ) +
  theme_minimal(base_family = font) +
  theme(
    legend.position = "none",
    axis.title.x = element_text(margin = margin(t = 10)),
    axis.text.x = element_text( size=15),
    panel.background = element_rect(fill = "white"),
    panel.grid.major = element_line(color = "gray80", size = 0.5),
    panel.grid.minor = element_line(color = "gray95", size = 0.3)
  ) +
  scale_y_continuous(expand = expansion(mult = c(0.1, 0.2)))


# Combine plots
combined_plot <- plot_proportion / plot_absolute +
  plot_layout(heights = c(1, 0.5)) &
  theme(plot.margin = margin(5,5, 0, 5))

# Display combined plot
combined_plot

2.3 PCA Plot

2.3.1 PladB

# Subset and scale data
pca_pladb <- pladb_data[, c("EVENT",pladb_events$Samples)] %>%
  na.omit()
pca_spire <- spire_data[,  c("EVENT",spire_events$Samples)]
pca_fmndko <- fmndko_data[,  c("EVENT",fmndko_events$Samples)]
rownames(pca_pladb)<- pca_pladb$EVENT
pca_pladb<-t(pca_pladb[,-1])


pca_result_pladb <- prcomp(pca_pladb)  # Perform PCA

# Prepare PCA results for plotting
pca_data <- as.data.frame(pca_result_pladb$x)
pca_data$Sample <- rownames(pca_data)
pca_data$condition<-metadata_pladb$Description

# Create PCA plot
pca_plot <- ggplot(pca_data, aes(x = PC1, y = PC2, color = as.factor(condition), )) +
  geom_point(size = 6) +
  labs(
    title = "<b style='font-size:18px;'>PCA of Samples (PSI Data)</b>",
    x = paste("PC1 (", round(100 * summary(pca_result_pladb)$importance[2, 1], 1), "%)", sep = ""),
    y = paste("PC2 (", round(100 * summary(pca_result_pladb)$importance[2, 2], 1), "%)", sep = ""),
    caption = paste0("Created by AG on ", Sys.Date()), size=3
  ) +
  theme_minimal(base_family = "roboto") +
  theme(
    plot.title = element_markdown(),
    plot.title.position = "plot",
    axis.title = element_text(size = 14),
    axis.text = element_text(size = 12),
    legend.position = "right",
    legend.title = element_blank(),
    panel.background = element_rect(fill = "gray98"),
    panel.grid = element_line(color = "gray90", size = 0.5),
    panel.grid.major = element_line(color = "gray90", size = 0.5),
    panel.grid.minor = element_blank()
  ) +
  scale_colour_paletteer_d("MetBrewer::Hokusai3") +
  coord_fixed()

pca_plot

2.3.2 SpireDKO

# Subset and scale data
pca_spire <- spire_data[, c("EVENT",spire_events$Samples)] %>%
  na.omit()

rownames(pca_spire)<- pca_spire$EVENT
pca_spire<-t(pca_spire[,-1])


pca_result_spire <- prcomp(pca_spire)  # Perform PCA

# Prepare PCA results for plotting
pca_data <- as.data.frame(pca_result_spire$x)
pca_data$Sample <- rownames(pca_data)
pca_data$condition<-metadata_spire$condition

# Create PCA plot
pca_plot <- ggplot(pca_data, aes(x = PC1, y = PC2, color = as.factor(condition), )) +
  geom_point(size = 6) +
  labs(
    title = "<b style='font-size:18px;'>PCA of Samples (PSI Data)</b>",
    x = paste("PC1 (", round(100 * summary(pca_result_pladb)$importance[2, 1], 1), "%)", sep = ""),
    y = paste("PC2 (", round(100 * summary(pca_result_pladb)$importance[2, 2], 1), "%)", sep = ""),
    caption = paste0("Created by AG on ", Sys.Date()), size=3
  ) +
  theme_minimal(base_family = "roboto") +
  theme(
    plot.title = element_markdown(),
    plot.title.position = "plot",
    axis.title = element_text(size = 14),
    axis.text = element_text(size = 12),
    legend.position = "right",
    legend.title = element_blank(),
    panel.background = element_rect(fill = "gray98"),
    panel.grid = element_line(color = "gray90", size = 0.5),
    panel.grid.major = element_line(color = "gray90", size = 0.5),
    panel.grid.minor = element_blank()
  ) +
  scale_colour_paletteer_d("MetBrewer::Hokusai3") +
  coord_fixed()

pca_plot

2.3.3 FMN2DKO

# Subset and scale data
pca_fmndko <- fmndko_data[, c("EVENT",fmndko_events$Samples)] %>%
  na.omit()

rownames(pca_fmndko)<- pca_fmndko$EVENT
pca_fmndko<-t(pca_fmndko[,-1])


pca_result_fmndko <- prcomp(pca_fmndko)  # Perform PCA

# Prepare PCA results for plotting
pca_data <- as.data.frame(pca_result_fmndko$x)
pca_data$Sample <- rownames(pca_data)
pca_data$condition<-metadata_fmn2dko$condition

# Create PCA plot
pca_plot <- ggplot(pca_data, aes(x = PC1, y = PC2, color = as.factor(condition), )) +
  geom_point(size = 6) +
  labs(
    title = "<b style='font-size:18px;'>PCA of Samples (PSI Data)</b>",
    x = paste("PC1 (", round(100 * summary(pca_result_pladb)$importance[2, 1], 1), "%)", sep = ""),
    y = paste("PC2 (", round(100 * summary(pca_result_pladb)$importance[2, 2], 1), "%)", sep = ""),
    caption = paste0("Created by AG on ", Sys.Date()), size=3
  ) +
  theme_minimal(base_family = "roboto") +
  theme(
    plot.title = element_markdown(),
    plot.title.position = "plot",
    axis.title = element_text(size = 14),
    axis.text = element_text(size = 12),
    legend.position = "right",
    legend.title = element_blank(),
    panel.background = element_rect(fill = "gray98"),
    panel.grid = element_line(color = "gray90", size = 0.5),
    panel.grid.major = element_line(color = "gray90", size = 0.5),
    panel.grid.minor = element_blank()
  ) +
  scale_colour_paletteer_d("MetBrewer::Hokusai3") +
  coord_fixed()

pca_plot

3 Differential Splicing Calculation

# Extract unique groups and sample IDs
# pladb

groupingVariable <- "Description"


groups_pladb <- unique(metadata_pladb[, groupingVariable])
samples_pladb <- metadata_pladb$fastq_files
samples_pladb <- paste0(
  "X",
  sub("\\.fastq\\.gz$", "", samples_pladb),
  "_merged_trimmed"
)

# fmn2dko

metadata_fmn2dko<-data.frame(samples=c("SRR6026682_SRR6026683_SRR6026684_merged_trimmed","SRR6026685_SRR6026686_SRR6026687_merged_trimmed","SRR6026688_SRR6026689_SRR6026690_merged_trimmed","SRR6026691_SRR6026692_SRR6026693_merged_trimmed"),condition=c("control","control","fmndko","fmndko"))
groups_fmndko <- unique(metadata_fmn2dko[, "condition"])
samples_fmndko <- metadata_fmn2dko$samples

# spiredko

groups_spire <- unique(metadata_spire[, "condition"])
samples_spire <- metadata_spire$samples


random_colors=c("#D8D97AFF", "#95C36EFF", "#74C8C3FF", "#5A97C1FF", "#295384FF", "#0A2E57FF")

# Create group list with metadata
groupList_pladb <- lapply(1:length(groups_pladb), function(i) {
  list(
    name = groups_pladb[i],
    samples = samples_pladb[metadata_pladb[, groupingVariable] == groups_pladb[i]],
    color = random_colors[i]
  )
})
names(groupList_pladb) <- groups_pladb

# Create group list with metadata
groupList_fmndko <- lapply(1:length(groups_fmndko), function(i) {
  list(
    name = groups_fmndko[i],
    samples = samples_fmndko[metadata_fmn2dko[, "condition"] == groups_fmndko[i]],
    color = random_colors[i]
  )
})
names(groupList_fmndko) <- groups_fmndko

# Create group list with metadata
groupList_spire <- lapply(1:length(groups_spire), function(i) {
  list(
    name = groups_spire[i],
    samples = samples_spire[metadata_spire[, "condition"] == groups_spire[i]],
    color = random_colors[i]
  )
})
names(groupList_spire) <- groups_spire

# Define groups
groupA_pladb    <- "ControlFmn2+-"
groupC_pladb    <- "Fmn2+-_+1mMPlatB"
groupB_pladb <- "Fmn2+-_+10mMPlatB"

groupA_fmndko<-"control"
groupB_fmndko<-"fmndko"

groupA_spire<-"control"
groupB_spire<-"spiredko"

# Define samples inside each group
samplesA_pladb    <- groupList_pladb[[groupA_pladb]]$samples
samplesB_pladb    <- groupList_pladb[[groupB_pladb]]$samples
samplesC_pladb    <- groupList_pladb[[groupC_pladb]]$samples
colsGroupA_pladb    <- convertCols(pladb_exons$PSI, samplesA_pladb)
colsGroupB_pladb   <- convertCols(pladb_exons$PSI, samplesB_pladb)
colsGroupC_pladb   <- convertCols(pladb_exons$PSI, samplesC_pladb)

samplesA_fmndko <- groupList_fmndko[[groupA_fmndko]]$samples
samplesB_fmndko <- groupList_fmndko[[groupB_fmndko]]$samples
colsGroupA_fmndko <- convertCols(fmndko_exons$PSI, samplesA_fmndko)
colsGroupB_fmndko <- convertCols(fmndko_exons$PSI, samplesB_fmndko)

samplesA_spire <- groupList_spire[[groupA_spire]]$samples
samplesB_spire <- groupList_spire[[groupB_spire]]$samples
colsGroupA_spire <- convertCols(spire_exons$PSI, samplesA_spire)
colsGroupB_spire <- convertCols(spire_exons$PSI, samplesB_spire)


set.seed(42) #Setting seed for downstream simulations of the beta distribution

3.1 Exon Lists

3.1.1 Calculations

pladb_pdiff_exons <- prepareTableVolcano(
  psitable = pladb_exons$PSI,
  qualtable = pladb_exons$Qual,
  npoints = 500,
  colsA = colsGroupA_pladb,
  colsB = colsGroupB_pladb,
  labA = groupA_pladb,
  labB = groupB_pladb,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed=TRUE,
  CoverageWeight = FALSE)

pladb_1mm_control_pdiff_exons <- prepareTableVolcano(
  psitable = pladb_exons$PSI,
  qualtable = pladb_exons$Qual,
  npoints = 500,
  colsA = colsGroupA_pladb,
  colsB = colsGroupC_pladb,
  labA = groupA_pladb,
  labB = groupC_pladb,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed=TRUE,
  CoverageWeight = FALSE)

pladb_1mm_10mm_pdiff_exons <- prepareTableVolcano(
  psitable = pladb_exons$PSI,
  qualtable = pladb_exons$Qual,
  npoints = 500,
  colsA = colsGroupC_pladb,
  colsB = colsGroupB_pladb,
  labA = groupC_pladb,
  labB = groupB_pladb,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed=TRUE,
  CoverageWeight = FALSE)

# Prepare table for fmn2dko
fmndko_pdiff_exons <- prepareTableVolcano(
  psitable = fmndko_exons$PSI,
  qualtable = fmndko_exons$Qual,
  npoints = 500,
  colsA = colsGroupA_fmndko,
  colsB = colsGroupB_fmndko,
  labA = groupA_fmndko,
  labB = groupB_fmndko,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed = TRUE,
  CoverageWeight = FALSE
)

# Prepare table for spire
spire_pdiff_exons <- prepareTableVolcano(
  psitable = spire_exons$PSI,
  qualtable = spire_exons$Qual,
  npoints = 500,
  colsA = colsGroupA_spire,
  colsB = colsGroupB_spire,
  labA = groupA_spire,
  labB = groupB_spire,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed = TRUE,
  CoverageWeight = FALSE
)

pladb_fstat_exons <- prepareTableVolcanoFstat(
  psitable = pladb_exons$PSI,
  qualtable = pladb_exons$Qual,
  npoints = 500,
  colsA = colsGroupA_pladb,
  colsB = colsGroupB_pladb,
  labA = groupA_pladb,
  labB = groupB_pladb,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed=TRUE,
  CoverageWeight = FALSE)

pladb_1mm_control_fstat_exons <- prepareTableVolcanoFstat(
  psitable = pladb_exons$PSI,
  qualtable = pladb_exons$Qual,
  npoints = 500,
  colsA = colsGroupA_pladb,
  colsB = colsGroupC_pladb,
  labA = groupA_pladb,
  labB = groupC_pladb,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed=TRUE,
  CoverageWeight = FALSE)

pladb_1mm_10mm_fstat_exons <- prepareTableVolcanoFstat(
  psitable = pladb_exons$PSI,
  qualtable = pladb_exons$Qual,
  npoints = 500,
  colsA = colsGroupC_pladb,
  colsB = colsGroupB_pladb,
  labA = groupC_pladb,
  labB = groupB_pladb,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed=TRUE,
  CoverageWeight = FALSE)

# Prepare table for fmndko
fmndko_fstat_exons <- prepareTableVolcanoFstat(
  psitable = fmndko_exons$PSI,
  qualtable = fmndko_exons$Qual,
  npoints = 500,
  colsA = colsGroupA_fmndko,
  colsB = colsGroupB_fmndko,
  labA = groupA_fmndko,
  labB = groupB_fmndko,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed = TRUE,
  CoverageWeight = FALSE
)


# Prepare table for spire
spire_fstat_exons <- prepareTableVolcanoFstat(
  psitable = spire_exons$PSI,
  qualtable = spire_exons$Qual,
  npoints = 500,
  colsA = colsGroupA_spire,
  colsB = colsGroupB_spire,
  labA = groupA_spire,
  labB = groupB_spire,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed = TRUE,
  CoverageWeight = FALSE
)

pladb_fdr_exons <- prepareTableVolcanoFDR(
  psitable = pladb_exons$PSI,
  qualtable = pladb_exons$Qual,
  npoints = 500,
  colsA = colsGroupA_pladb,
  colsB = colsGroupB_pladb,
  labA = groupA_pladb,
  labB = groupB_pladb,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  nsim = 100,
  seed=TRUE,
  CoverageWeight = FALSE)

pladb_1mm_control_fdr_exons <- prepareTableVolcanoFDR(
  psitable = pladb_exons$PSI,
  qualtable = pladb_exons$Qual,
  npoints = 500,
  colsA = colsGroupA_pladb,
  colsB = colsGroupC_pladb,
  labA = groupA_pladb,
  labB = groupC_pladb,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  nsim = 100,
  seed=TRUE,
  CoverageWeight = FALSE)


pladb_1mm_10mm_fdr_exons <- prepareTableVolcanoFDR(
  psitable = pladb_exons$PSI,
  qualtable = pladb_exons$Qual,
  npoints = 500,
  colsA = colsGroupC_pladb,
  colsB = colsGroupB_pladb,
  labA = groupC_pladb,
  labB = groupB_pladb,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  nsim = 100,
  seed=TRUE,
  CoverageWeight = FALSE)


# Prepare table for fmndko
fmndko_fdr_exons <- prepareTableVolcanoFDR(
  psitable = fmndko_exons$PSI,
  qualtable = fmndko_exons$Qual,
  npoints = 500,
  colsA = colsGroupA_fmndko,
  colsB = colsGroupB_fmndko,
  labA = groupA_fmndko,
  labB = groupB_fmndko,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  nsim = 100,
  seed = TRUE,
  CoverageWeight = FALSE
)

# Prepare table for spire
spire_fdr_exons <- prepareTableVolcanoFDR(
  psitable = spire_exons$PSI,
  qualtable = spire_exons$Qual,
  npoints = 500,
  colsA = colsGroupA_spire,
  colsB = colsGroupB_spire,
  labA = groupA_spire,
  labB = groupB_spire,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  nsim = 100,
  seed = TRUE,
  CoverageWeight = FALSE
)

3.1.2 FDR Volcano Plots

3.1.2.1 PladB 10mM vs Control

plotVolcanoFDR(betasTable =filter(pladb_fdr_exons,!is.na(EVENT)),
                            labA = groupA_pladb,
                            labB = groupB_pladb,
                            basalColor = "#89C0AE",
                            interestColor = "#E69A9C") +

theme(
  plot.title = element_text(size = 10),  # Title font size
  axis.title = element_text(size = 10),  # Axis title font size
  axis.text = element_text(size = 10),   # Axis text font size
  legend.text = element_text(size = 10), # Legend text font size
  legend.title = element_text(size = 10) # Legend title font size
)

3.1.2.2 PladB 1mM vs Control

plotVolcanoFDR(betasTable =filter(pladb_1mm_control_fdr_exons,!is.na(EVENT)),
                            labA = groupA_pladb,
                            labB = groupC_pladb,
                            basalColor = "#89C0AE",
                            interestColor = "#E69A9C") +

theme(
  plot.title = element_text(size = 10),  # Title font size
  axis.title = element_text(size = 10),  # Axis title font size
  axis.text = element_text(size = 10),   # Axis text font size
  legend.text = element_text(size = 10), # Legend text font size
  legend.title = element_text(size = 10) # Legend title font size
)

3.1.2.3 PladB 10mM vs 1mM

plotVolcanoFDR(betasTable =filter(pladb_1mm_10mm_fdr_exons,!is.na(EVENT)),
                            labA = groupC_pladb,
                            labB = groupB_pladb,
                            basalColor = "#89C0AE",
                            interestColor = "#E69A9C") +

theme(
  plot.title = element_text(size = 10),  # Title font size
  axis.title = element_text(size = 10),  # Axis title font size
  axis.text = element_text(size = 10),   # Axis text font size
  legend.text = element_text(size = 10), # Legend text font size
  legend.title = element_text(size = 10) # Legend title font size
)

3.1.2.4 FMN2DKO

plotVolcanoFDR(betasTable =filter(fmndko_fdr_exons,!is.na(EVENT)),
                            labA = groupA_fmndko,
                            labB = groupB_fmndko,
                            basalColor = "#89C0AE",
                            interestColor = "#E69A9C") +
theme(
  plot.title = element_text(size = 10),  # Title font size
  axis.title = element_text(size = 10),  # Axis title font size
  axis.text = element_text(size = 10),   # Axis text font size
  legend.text = element_text(size = 10), # Legend text font size
  legend.title = element_text(size = 10) # Legend title font size
)

3.1.2.5 Spire2

plotVolcanoFDR(betasTable =filter(spire_fdr_exons,!is.na(EVENT)),
                            labA = groupA_spire,
                            labB = groupB_spire,
                            basalColor = "#89C0AE",
                            interestColor = "#E69A9C") +
theme(
  plot.title = element_text(size = 10),  # Title font size
  axis.title = element_text(size = 10),  # Axis title font size
  axis.text = element_text(size = 10),   # Axis text font size
  legend.text = element_text(size = 10), # Legend text font size
  legend.title = element_text(size = 10) # Legend title font size
)

3.1.3 Individual Combination Tables

3.1.3.1 PladB 10mM vs control

pladb_pdiff_exons_tab <- pladb_pdiff_exons[, c("GENE", "EVENT", "COORD", "Pdiff", "deltapsi")] %>%
  filter(!is.na(EVENT), Pdiff >= 0.95) %>%  # Remove rows where EVENT is NA and Pdiff >= 0.95
  arrange(desc(Pdiff), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(pladb_pdiff_exons_tab) <- 1:nrow(pladb_pdiff_exons_tab)

pladb_fstat_exons_tab <- pladb_fstat_exons[, c("GENE","EVENT","COORD","Fstat","deltapsi")] %>%
  filter(!is.na(EVENT)) %>%  # Remove rows where EVENT is NA
  arrange(desc(Fstat), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(pladb_fstat_exons_tab) <- 1:nrow(pladb_fstat_exons_tab)

pladb_fdr_exons_tab <- pladb_fdr_exons[, c("GENE","EVENT","COORD","FDR","deltapsi")] %>%
  filter(!is.na(EVENT), FDR<=0.05) %>%  # Remove rows where EVENT is NA
  arrange((FDR), desc(abs(deltapsi)))
rownames(pladb_fdr_exons_tab) <- 1:nrow(pladb_fdr_exons_tab)

# Find common EVENTs across all tables
common_events <- Reduce(intersect, list(pladb_pdiff_exons_tab$EVENT, pladb_fstat_exons_tab$EVENT, pladb_fdr_exons_tab$EVENT))

# Filter each table to include only rows with common EVENTs
filtered_exons_diff_tab <- pladb_pdiff_exons_tab %>% filter(EVENT %in% common_events)
filtered_exons_fstat_tab <- pladb_fstat_exons_tab %>% filter(EVENT %in% common_events)
filtered_exons_fdr_tab <- pladb_fdr_exons_tab %>% filter(EVENT %in% common_events)


# Select only unique columns from each table
filtered_exons_diff_tab <- filtered_exons_diff_tab %>% select(GENE,EVENT,COORD, deltapsi,Pdiff) # Extra column is Pdiff
filtered_exons_fstat_tab <- filtered_exons_fstat_tab %>% select(EVENT, Fstat) # Extra column is Fstat
filtered_exons_fdr_tab <- filtered_exons_fdr_tab %>% select(EVENT, FDR) # Extra column is FDR

# Merge the data frames by EVENT
pladb_exon_combination_tab <- filtered_exons_diff_tab %>%
  inner_join(filtered_exons_fstat_tab, by = "EVENT") %>%
  inner_join(filtered_exons_fdr_tab, by = "EVENT") %>%
  distinct() %>%
  arrange(desc(abs(deltapsi)), desc(Fstat),desc(FDR), desc(Pdiff))

# Render DataTable with enhancements
datatable(
  pladb_exon_combination_tab,
  options = list(
    pageLength = 10,                      # Rows per page
    autoWidth = TRUE,                     # Adjust column widths automatically
    dom = 'Bfrtip',                       # Add buttons for export
    buttons = c("copy", "csv", "excel", "pdf"),  # Simplified button definitions
    columnDefs = list(
      list(targets = "_all", className = "dt-center"),    # Center-align all columns
      list(
        targets = 3,  # Highlight significant Pdiff
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && parseFloat(data) > 0.99) {
               return '<span style=\"color: green; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      ),
      list(
        targets = 4,  # Highlight large deltapsi values
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && Math.abs(parseFloat(data)) > 0.5) {
               return '<span style=\"color: red; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      )
    )
  ),
  rownames = FALSE,                       # Disable row names
  extensions = "Buttons"                  # Enable export buttons
) %>%
  formatStyle(
    columns = "Pdiff",
    color = styleInterval(0.99, c("black", "green")),
    fontWeight = styleEqual(0.99, "bold")  # Ensure consistent styling
  ) %>%
  formatStyle(
  columns = "deltapsi",
  color = styleInterval(c(-0.5, 0.5), c("black", "red", "black")),  # Three colors for two intervals
  fontWeight = styleInterval(c(-0.5, 0.5), c("normal", "bold", "normal"))  # Three weights for two intervals
)

write.csv(pladb_exon_combination_tab,
          file = paste0("pladb_exon_combination_table_", format(Sys.Date(), "%Y%m%d"), ".csv"),
          row.names = FALSE)

3.1.3.2 PladB 1mM vs Control

pladb_pdiff_exons_tab <- pladb_1mm_control_pdiff_exons[, c("GENE", "EVENT", "COORD", "Pdiff", "deltapsi")] %>%
  filter(!is.na(EVENT), Pdiff >= 0.95) %>%  # Remove rows where EVENT is NA and Pdiff >= 0.95
  arrange(desc(Pdiff), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(pladb_pdiff_exons_tab) <- 1:nrow(pladb_pdiff_exons_tab)

pladb_fstat_exons_tab <- pladb_1mm_control_fstat_exons[, c("GENE","EVENT","COORD","Fstat","deltapsi")] %>%
  filter(!is.na(EVENT)) %>%  # Remove rows where EVENT is NA
  arrange(desc(Fstat), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(pladb_fstat_exons_tab) <- 1:nrow(pladb_fstat_exons_tab)

pladb_fdr_exons_tab <- pladb_1mm_control_fdr_exons[, c("GENE","EVENT","COORD","FDR","deltapsi")] %>%
  filter(!is.na(EVENT), FDR<=0.05) %>%  # Remove rows where EVENT is NA
  arrange((FDR), desc(abs(deltapsi)))
rownames(pladb_fdr_exons_tab) <- 1:nrow(pladb_fdr_exons_tab)

# Find common EVENTs across all tables
common_events <- Reduce(intersect, list(pladb_pdiff_exons_tab$EVENT, pladb_fstat_exons_tab$EVENT, pladb_fdr_exons_tab$EVENT))

# Filter each table to include only rows with common EVENTs
filtered_exons_diff_tab <- pladb_pdiff_exons_tab %>% filter(EVENT %in% common_events)
filtered_exons_fstat_tab <- pladb_fstat_exons_tab %>% filter(EVENT %in% common_events)
filtered_exons_fdr_tab <- pladb_fdr_exons_tab %>% filter(EVENT %in% common_events)


# Select only unique columns from each table
filtered_exons_diff_tab <- filtered_exons_diff_tab %>% select(GENE,EVENT,COORD, deltapsi,Pdiff) # Extra column is Pdiff
filtered_exons_fstat_tab <- filtered_exons_fstat_tab %>% select(EVENT, Fstat) # Extra column is Fstat
filtered_exons_fdr_tab <- filtered_exons_fdr_tab %>% select(EVENT, FDR) # Extra column is FDR

# Merge the data frames by EVENT
pladb_1mm_control_exon_combination_tab <- filtered_exons_diff_tab %>%
  inner_join(filtered_exons_fstat_tab, by = "EVENT") %>%
  inner_join(filtered_exons_fdr_tab, by = "EVENT") %>%
  distinct() %>%
  arrange(desc(abs(deltapsi)), desc(Fstat),desc(FDR), desc(Pdiff))

# Find if events appear in the control to 10mM pairwise comparison
pladb_1mm_control_exon_combination_tab$appears_in_10mm_control <- pladb_1mm_control_exon_combination_tab$EVENT %in% pladb_exon_combination_tab$EVENT

# Render DataTable with enhancements
datatable(
  pladb_1mm_control_exon_combination_tab,
  options = list(
    pageLength = 10,                      # Rows per page
    autoWidth = TRUE,                     # Adjust column widths automatically
    dom = 'Bfrtip',                       # Add buttons for export
    buttons = c("copy", "csv", "excel", "pdf"),  # Simplified button definitions
    columnDefs = list(
      list(targets = "_all", className = "dt-center"),    # Center-align all columns
      list(
        targets = 3,  # Highlight significant Pdiff
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && parseFloat(data) > 0.99) {
               return '<span style=\"color: green; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      ),
      list(
        targets = 4,  # Highlight large deltapsi values
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && Math.abs(parseFloat(data)) > 0.5) {
               return '<span style=\"color: red; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      )
    )
  ),
  rownames = FALSE,                       # Disable row names
  extensions = "Buttons"                  # Enable export buttons
) %>%
  formatStyle(
    columns = "Pdiff",
    color = styleInterval(0.99, c("black", "green")),
    fontWeight = styleEqual(0.99, "bold")  # Ensure consistent styling
  ) %>%
  formatStyle(
  columns = "deltapsi",
  color = styleInterval(c(-0.5, 0.5), c("black", "red", "black")),  # Three colors for two intervals
  fontWeight = styleInterval(c(-0.5, 0.5), c("normal", "bold", "normal"))  # Three weights for two intervals
)

write.csv(pladb_1mm_control_exon_combination_tab,
          file = paste0("pladb_1mm_control_exon_combination_table_", format(Sys.Date(), "%Y%m%d"), ".csv"),
          row.names = FALSE)

3.1.3.3 PladB 10mM vs 1mM

pladb_pdiff_exons_tab <- pladb_1mm_10mm_pdiff_exons[, c("GENE", "EVENT", "COORD", "Pdiff", "deltapsi")] %>%
  filter(!is.na(EVENT), Pdiff >= 0.95) %>%  # Remove rows where EVENT is NA and Pdiff >= 0.95
  arrange(desc(Pdiff), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(pladb_pdiff_exons_tab) <- 1:nrow(pladb_pdiff_exons_tab)

pladb_fstat_exons_tab <- pladb_1mm_10mm_fstat_exons[, c("GENE","EVENT","COORD","Fstat","deltapsi")] %>%
  filter(!is.na(EVENT)) %>%  # Remove rows where EVENT is NA
  arrange(desc(Fstat), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(pladb_fstat_exons_tab) <- 1:nrow(pladb_fstat_exons_tab)

pladb_fdr_exons_tab <- pladb_1mm_10mm_fdr_exons[, c("GENE","EVENT","COORD","FDR","deltapsi")] %>%
  filter(!is.na(EVENT), FDR<=0.05) %>%  # Remove rows where EVENT is NA
  arrange((FDR), desc(abs(deltapsi)))
rownames(pladb_fdr_exons_tab) <- 1:nrow(pladb_fdr_exons_tab)

# Find common EVENTs across all tables
common_events <- Reduce(intersect, list(pladb_pdiff_exons_tab$EVENT, pladb_fstat_exons_tab$EVENT, pladb_fdr_exons_tab$EVENT))

# Filter each table to include only rows with common EVENTs
filtered_exons_diff_tab <- pladb_pdiff_exons_tab %>% filter(EVENT %in% common_events)
filtered_exons_fstat_tab <- pladb_fstat_exons_tab %>% filter(EVENT %in% common_events)
filtered_exons_fdr_tab <- pladb_fdr_exons_tab %>% filter(EVENT %in% common_events)


# Select only unique columns from each table
filtered_exons_diff_tab <- filtered_exons_diff_tab %>% select(GENE,EVENT,COORD, deltapsi,Pdiff) # Extra column is Pdiff
filtered_exons_fstat_tab <- filtered_exons_fstat_tab %>% select(EVENT, Fstat) # Extra column is Fstat
filtered_exons_fdr_tab <- filtered_exons_fdr_tab %>% select(EVENT, FDR) # Extra column is FDR

# Merge the data frames by EVENT
pladb_10mm_1mm_exon_combination_tab <- filtered_exons_diff_tab %>%
  inner_join(filtered_exons_fstat_tab, by = "EVENT") %>%
  inner_join(filtered_exons_fdr_tab, by = "EVENT") %>%
  distinct() %>%
  arrange(desc(abs(deltapsi)), desc(Fstat),desc(FDR), desc(Pdiff))

# Find if events appear in the control to 10mM pairwise comparison
pladb_10mm_1mm_exon_combination_tab$appears_in_10mm_control <- pladb_10mm_1mm_exon_combination_tab$EVENT %in% pladb_exon_combination_tab$EVENT

# Render DataTable with enhancements
datatable(
  pladb_10mm_1mm_exon_combination_tab,
  options = list(
    pageLength = 10,                      # Rows per page
    autoWidth = TRUE,                     # Adjust column widths automatically
    dom = 'Bfrtip',                       # Add buttons for export
    buttons = c("copy", "csv", "excel", "pdf"),  # Simplified button definitions
    columnDefs = list(
      list(targets = "_all", className = "dt-center"),    # Center-align all columns
      list(
        targets = 3,  # Highlight significant Pdiff
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && parseFloat(data) > 0.99) {
               return '<span style=\"color: green; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      ),
      list(
        targets = 4,  # Highlight large deltapsi values
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && Math.abs(parseFloat(data)) > 0.5) {
               return '<span style=\"color: red; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      )
    )
  ),
  rownames = FALSE,                       # Disable row names
  extensions = "Buttons"                  # Enable export buttons
) %>%
  formatStyle(
    columns = "Pdiff",
    color = styleInterval(0.99, c("black", "green")),
    fontWeight = styleEqual(0.99, "bold")  # Ensure consistent styling
  ) %>%
  formatStyle(
  columns = "deltapsi",
  color = styleInterval(c(-0.5, 0.5), c("black", "red", "black")),  # Three colors for two intervals
  fontWeight = styleInterval(c(-0.5, 0.5), c("normal", "bold", "normal"))  # Three weights for two intervals
)

write.csv(pladb_10mm_1mm_exon_combination_tab,
          file = paste0("pladb_10mm_1mm_exon_combination_table_", format(Sys.Date(), "%Y%m%d"), ".csv"),
          row.names = FALSE)

3.1.3.4 FMN2DKO

fmndko_pdiff_exons_tab <- fmndko_pdiff_exons[, c("GENE", "EVENT", "COORD", "Pdiff", "deltapsi")] %>%
  filter(!is.na(EVENT), Pdiff >= 0.95) %>%  # Remove rows where EVENT is NA and Pdiff >= 0.95
  arrange(desc(Pdiff), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(fmndko_pdiff_exons_tab) <- 1:nrow(fmndko_pdiff_exons_tab)

fmndko_fstat_exons_tab <- fmndko_fstat_exons[, c("GENE","EVENT","COORD","Fstat","deltapsi")] %>%
  filter(!is.na(EVENT)) %>%  # Remove rows where EVENT is NA
  arrange(desc(Fstat), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(fmndko_fstat_exons_tab) <- 1:nrow(fmndko_fstat_exons_tab)

fmndko_fdr_exons_tab <- fmndko_fdr_exons[, c("GENE","EVENT","COORD","FDR","deltapsi")] %>%
  filter(!is.na(EVENT), FDR<=0.05) %>%  # Remove rows where EVENT is NA
  arrange((FDR), desc(abs(deltapsi)))
rownames(fmndko_fdr_exons_tab) <- 1:nrow(fmndko_fdr_exons_tab)

# Find common EVENTs across all tables
common_events <- Reduce(intersect, list(fmndko_pdiff_exons_tab$EVENT, fmndko_fstat_exons_tab$EVENT, fmndko_fdr_exons_tab$EVENT))

# Filter each table to include only rows with common EVENTs
filtered_exon_diff_tab <- fmndko_pdiff_exons_tab %>% filter(EVENT %in% common_events)
filtered_exon_fstat_tab <- fmndko_fstat_exons_tab %>% filter(EVENT %in% common_events)
filtered_exon_fdr_tab <- fmndko_fdr_exons_tab %>% filter(EVENT %in% common_events)


# Select only unique columns from each table
filtered_exon_diff_tab <- filtered_exon_diff_tab %>% select(GENE,EVENT,COORD, deltapsi,Pdiff) # Extra column is Pdiff
filtered_exon_fstat_tab <- filtered_exon_fstat_tab %>% select(EVENT, Fstat) # Extra column is Fstat
filtered_exon_fdr_tab <- filtered_exon_fdr_tab %>% select(EVENT, FDR) # Extra column is FDR

# Merge the data frames by EVENT
fmndko_exon_combination_tab <- filtered_exon_diff_tab %>%
  inner_join(filtered_exon_fstat_tab, by = "EVENT") %>%
  inner_join(filtered_exon_fdr_tab, by = "EVENT") %>%
  distinct() %>%
  arrange(desc(abs(deltapsi)), desc(Fstat),desc(FDR), desc(Pdiff))

# Render DataTable with enhancements
datatable(
  fmndko_exon_combination_tab,
  options = list(
    pageLength = 10,                      # Rows per page
    autoWidth = TRUE,                     # Adjust column widths automatically
    dom = 'Bfrtip',                       # Add buttons for export
    buttons = c("copy", "csv", "excel", "pdf"),  # Simplified button definitions
    columnDefs = list(
      list(targets = "_all", className = "dt-center"),    # Center-align all columns
      list(
        targets = 3,  # Highlight significant Pdiff
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && parseFloat(data) > 0.99) {
               return '<span style=\"color: green; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      ),
      list(
        targets = 4,  # Highlight large deltapsi values
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && Math.abs(parseFloat(data)) > 0.5) {
               return '<span style=\"color: red; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      )
    )
  ),
  rownames = FALSE,                       # Disable row names
  extensions = "Buttons"                  # Enable export buttons
) %>%
  formatStyle(
    columns = "Pdiff",
    color = styleInterval(0.99, c("black", "green")),
    fontWeight = styleEqual(0.99, "bold")  # Ensure consistent styling
  ) %>%
  formatStyle(
  columns = "deltapsi",
  color = styleInterval(c(-0.5, 0.5), c("black", "red", "black")),  # Three colors for two intervals
  fontWeight = styleInterval(c(-0.5, 0.5), c("normal", "bold", "normal"))  # Three weights for two intervals
)

write.csv(fmndko_exon_combination_tab,
          file = paste0("fmndko_exon_combination_table_", format(Sys.Date(), "%Y%m%d"), ".csv"),
          row.names = FALSE)

3.1.3.5 Spire2

spire_pdiff_exons_tab <- spire_pdiff_exons[, c("GENE", "EVENT", "COORD", "Pdiff", "deltapsi")] %>%
  filter(!is.na(EVENT), Pdiff >= 0.95) %>%  # Remove rows where EVENT is NA and Pdiff >= 0.95
  arrange(desc(Pdiff), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(spire_pdiff_exons_tab) <- 1:nrow(spire_pdiff_exons_tab)

spire_fstat_exons_tab <- spire_fstat_exons[, c("GENE","EVENT","COORD","Fstat","deltapsi")] %>%
  filter(!is.na(EVENT)) %>%  # Remove rows where EVENT is NA
  arrange(desc(Fstat), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(spire_fstat_exons_tab) <- 1:nrow(spire_fstat_exons_tab)

spire_fdr_exons_tab <- spire_fdr_exons[, c("GENE","EVENT","COORD","FDR","deltapsi")] %>%
  filter(!is.na(EVENT), FDR<=0.05) %>%  # Remove rows where EVENT is NA
  arrange((FDR), desc(abs(deltapsi)))
rownames(spire_fdr_exons_tab) <- 1:nrow(spire_fdr_exons_tab)

# Find common EVENTs across all tables
common_events <- Reduce(intersect, list(spire_pdiff_exons_tab$EVENT, spire_fstat_exons_tab$EVENT, spire_fdr_exons_tab$EVENT))

# Filter each table to include only rows with common EVENTs
filtered_exon_diff_tab <- spire_pdiff_exons_tab %>% filter(EVENT %in% common_events)
filtered_exon_fstat_tab <- spire_fstat_exons_tab %>% filter(EVENT %in% common_events)
filtered_exon_fdr_tab <- spire_fdr_exons_tab %>% filter(EVENT %in% common_events)


# Select only unique columns from each table
filtered_exon_diff_tab <- filtered_exon_diff_tab %>% select(GENE,EVENT,COORD, deltapsi,Pdiff) # Extra column is Pdiff
filtered_exon_fstat_tab <- filtered_exon_fstat_tab %>% select(EVENT, Fstat) # Extra column is Fstat
filtered_exon_fdr_tab <- filtered_exon_fdr_tab %>% select(EVENT, FDR) # Extra column is FDR

# Merge the data frames by EVENT
spire_exon_combination_tab <- filtered_exon_diff_tab %>%
  inner_join(filtered_exon_fstat_tab, by = "EVENT") %>%
  inner_join(filtered_exon_fdr_tab, by = "EVENT") %>%
  distinct() %>%
  arrange(desc(abs(deltapsi)), desc(Fstat),desc(FDR), desc(Pdiff))

# Render DataTable with enhancements
datatable(
  spire_exon_combination_tab,
  options = list(
    pageLength = 10,                      # Rows per page
    autoWidth = TRUE,                     # Adjust column widths automatically
    dom = 'Bfrtip',                       # Add buttons for export
    buttons = c("copy", "csv", "excel", "pdf"),  # Simplified button definitions
    columnDefs = list(
      list(targets = "_all", className = "dt-center"),    # Center-align all columns
      list(
        targets = 3,  # Highlight significant Pdiff
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && parseFloat(data) > 0.99) {
               return '<span style=\"color: green; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      ),
      list(
        targets = 4,  # Highlight large deltapsi values
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && Math.abs(parseFloat(data)) > 0.5) {
               return '<span style=\"color: red; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      )
    )
  ),
  rownames = FALSE,                       # Disable row names
  extensions = "Buttons"                  # Enable export buttons
) %>%
  formatStyle(
    columns = "Pdiff",
    color = styleInterval(0.99, c("black", "green")),
    fontWeight = styleEqual(0.99, "bold")  # Ensure consistent styling
  ) %>%
  formatStyle(
  columns = "deltapsi",
  color = styleInterval(c(-0.5, 0.5), c("black", "red", "black")),  # Three colors for two intervals
  fontWeight = styleInterval(c(-0.5, 0.5), c("normal", "bold", "normal"))  # Three weights for two intervals
)

write.csv(spire_exon_combination_tab,
          file = paste0("spire_exon_combination_table_", format(Sys.Date(), "%Y%m%d"), ".csv"),
          row.names = FALSE)

3.2 Intron Lists

3.2.1 Calculations

pladb_pdiff_introns <- prepareTableVolcano(
  psitable = pladb_introns$PSI,
  qualtable = pladb_introns$Qual,
  npoints = 500,
  colsA = colsGroupA_pladb,
  colsB = colsGroupB_pladb,
  labA = groupA_pladb,
  labB = groupB_pladb,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed=TRUE,
  CoverageWeight = FALSE)

pladb_1mm_control_pdiff_introns <- prepareTableVolcano(
  psitable = pladb_introns$PSI,
  qualtable = pladb_introns$Qual,
  npoints = 500,
  colsA = colsGroupA_pladb,
  colsB = colsGroupC_pladb,
  labA = groupA_pladb,
  labB = groupC_pladb,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed=TRUE,
  CoverageWeight = FALSE)

pladb_1mm_10mm_pdiff_introns <- prepareTableVolcano(
  psitable = pladb_introns$PSI,
  qualtable = pladb_introns$Qual,
  npoints = 500,
  colsA = colsGroupC_pladb,
  colsB = colsGroupB_pladb,
  labA = groupC_pladb,
  labB = groupB_pladb,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed=TRUE,
  CoverageWeight = FALSE)

# Prepare table for Tao
fmndko_pdiff_introns <- prepareTableVolcano(
  psitable = fmndko_introns$PSI,
  qualtable = fmndko_introns$Qual,
  npoints = 500,
  colsA = colsGroupA_fmndko,
  colsB = colsGroupB_fmndko,
  labA = groupA_fmndko,
  labB = groupB_fmndko,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed = TRUE,
  CoverageWeight = FALSE
)

# Prepare table for spire
spire_pdiff_introns <- prepareTableVolcano(
  psitable = spire_introns$PSI,
  qualtable = spire_introns$Qual,
  npoints = 500,
  colsA = colsGroupA_spire,
  colsB = colsGroupB_spire,
  labA = groupA_spire,
  labB = groupB_spire,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed = TRUE,
  CoverageWeight = FALSE
)

pladb_fstat_introns <- prepareTableVolcanoFstat(
  psitable = pladb_introns$PSI,
  qualtable = pladb_introns$Qual,
  npoints = 500,
  colsA = colsGroupA_pladb,
  colsB = colsGroupB_pladb,
  labA = groupA_pladb,
  labB = groupB_pladb,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed=TRUE,
  CoverageWeight = FALSE)

pladb_1mm_control_fstat_introns <- prepareTableVolcanoFstat(
  psitable = pladb_introns$PSI,
  qualtable = pladb_introns$Qual,
  npoints = 500,
  colsA = colsGroupA_pladb,
  colsB = colsGroupC_pladb,
  labA = groupA_pladb,
  labB = groupC_pladb,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed=TRUE,
  CoverageWeight = FALSE)

pladb_1mm_10mm_fstat_introns <- prepareTableVolcanoFstat(
  psitable = pladb_introns$PSI,
  qualtable = pladb_introns$Qual,
  npoints = 500,
  colsA = colsGroupC_pladb,
  colsB = colsGroupB_pladb,
  labA = groupC_pladb,
  labB = groupB_pladb,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed=TRUE,
  CoverageWeight = FALSE)

# Prepare table for fmndko
fmndko_fstat_introns <- prepareTableVolcanoFstat(
  psitable = fmndko_introns$PSI,
  qualtable = fmndko_introns$Qual,
  npoints = 500,
  colsA = colsGroupA_fmndko,
  colsB = colsGroupB_fmndko,
  labA = groupA_fmndko,
  labB = groupB_fmndko,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed = TRUE,
  CoverageWeight = FALSE
)

# Prepare table for spire
spire_fstat_introns <- prepareTableVolcanoFstat(
  psitable = spire_introns$PSI,
  qualtable = spire_introns$Qual,
  npoints = 500,
  colsA = colsGroupA_spire,
  colsB = colsGroupB_spire,
  labA = groupA_spire,
  labB = groupB_spire,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed = TRUE,
  CoverageWeight = FALSE
)

pladb_fdr_introns <- prepareTableVolcanoFDR(
  psitable = pladb_introns$PSI,
  qualtable = pladb_introns$Qual,
  npoints = 500,
  colsA = colsGroupA_pladb,
  colsB = colsGroupB_pladb,
  labA = groupA_pladb,
  labB = groupB_pladb,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  nsim = 100,
  seed=TRUE,
  CoverageWeight = FALSE)

pladb_1mm_control_fdr_introns <- prepareTableVolcanoFDR(
  psitable = pladb_introns$PSI,
  qualtable = pladb_introns$Qual,
  npoints = 500,
  colsA = colsGroupA_pladb,
  colsB = colsGroupC_pladb,
  labA = groupA_pladb,
  labB = groupC_pladb,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  nsim = 100,
  seed=TRUE,
  CoverageWeight = FALSE)


pladb_1mm_10mm_fdr_introns <- prepareTableVolcanoFDR(
  psitable = pladb_introns$PSI,
  qualtable = pladb_introns$Qual,
  npoints = 500,
  colsA = colsGroupC_pladb,
  colsB = colsGroupB_pladb,
  labA = groupC_pladb,
  labB = groupB_pladb,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  nsim = 100,
  seed=TRUE,
  CoverageWeight = FALSE)


# Prepare table for fmndko
fmndko_fdr_introns <- prepareTableVolcanoFDR(
  psitable = fmndko_introns$PSI,
  qualtable = fmndko_introns$Qual,
  npoints = 500,
  colsA = colsGroupA_fmndko,
  colsB = colsGroupB_fmndko,
  labA = groupA_fmndko,
  labB = groupB_fmndko,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  nsim = 100,
  seed = TRUE,
  CoverageWeight = FALSE
)

# Prepare table for spire
spire_fdr_introns <- prepareTableVolcanoFDR(
  psitable = spire_introns$PSI,
  qualtable = spire_introns$Qual,
  npoints = 500,
  colsA = colsGroupA_spire,
  colsB = colsGroupB_spire,
  labA = groupA_spire,
  labB = groupB_spire,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  nsim = 100,
  seed = TRUE,
  CoverageWeight = FALSE
)

3.2.2 FDR Volcano Plots

3.2.2.1 PladB 10mM vs Control

plotVolcanoFDR(betasTable =filter(pladb_fdr_introns,!is.na(EVENT)),
                            labA = groupA_pladb,
                            labB = groupB_pladb,
                            basalColor = "#89C0AE",
                            interestColor = "#E69A9C") +

theme(
  plot.title = element_text(size = 10),  # Title font size
  axis.title = element_text(size = 10),  # Axis title font size
  axis.text = element_text(size = 10),   # Axis text font size
  legend.text = element_text(size = 10), # Legend text font size
  legend.title = element_text(size = 10) # Legend title font size
)

3.2.2.2 PladB 1mM vs Control

plotVolcanoFDR(betasTable =filter(pladb_1mm_control_fdr_introns,!is.na(EVENT)),
                            labA = groupA_pladb,
                            labB = groupC_pladb,
                            basalColor = "#89C0AE",
                            interestColor = "#E69A9C") +

theme(
  plot.title = element_text(size = 10),  # Title font size
  axis.title = element_text(size = 10),  # Axis title font size
  axis.text = element_text(size = 10),   # Axis text font size
  legend.text = element_text(size = 10), # Legend text font size
  legend.title = element_text(size = 10) # Legend title font size
)

3.2.2.3 PladB 10mM vs 1mM

plotVolcanoFDR(betasTable =filter(pladb_1mm_10mm_fdr_introns,!is.na(EVENT)),
                            labA = groupC_pladb,
                            labB = groupB_pladb,
                            basalColor = "#89C0AE",
                            interestColor = "#E69A9C") +

theme(
  plot.title = element_text(size = 10),  # Title font size
  axis.title = element_text(size = 10),  # Axis title font size
  axis.text = element_text(size = 10),   # Axis text font size
  legend.text = element_text(size = 10), # Legend text font size
  legend.title = element_text(size = 10) # Legend title font size
)

3.2.2.4 FMN2DKO

plotVolcanoFDR(betasTable =filter(fmndko_fdr_introns,!is.na(EVENT)),
                            labA = groupA_fmndko,
                            labB = groupB_fmndko,
                            basalColor = "#89C0AE",
                            interestColor = "#E69A9C") +
theme(
  plot.title = element_text(size = 10),  # Title font size
  axis.title = element_text(size = 10),  # Axis title font size
  axis.text = element_text(size = 10),   # Axis text font size
  legend.text = element_text(size = 10), # Legend text font size
  legend.title = element_text(size = 10) # Legend title font size
)

3.2.2.5 Spire2

plotVolcanoFDR(betasTable =filter(spire_fdr_introns,!is.na(EVENT)),
                            labA = groupA_spire,
                            labB = groupB_spire,
                            basalColor = "#89C0AE",
                            interestColor = "#E69A9C") +
theme(
  plot.title = element_text(size = 10),  # Title font size
  axis.title = element_text(size = 10),  # Axis title font size
  axis.text = element_text(size = 10),   # Axis text font size
  legend.text = element_text(size = 10), # Legend text font size
  legend.title = element_text(size = 10) # Legend title font size
)

3.2.3 Individual Combination Tables

3.2.3.1 PladB 10mM vs control

pladb_pdiff_introns_tab <- pladb_pdiff_introns[, c("GENE", "EVENT", "COORD", "Pdiff", "deltapsi")] %>%
  filter(!is.na(EVENT), Pdiff >= 0.95) %>%  # Remove rows where EVENT is NA and Pdiff >= 0.95
  arrange(desc(Pdiff), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(pladb_pdiff_introns_tab) <- 1:nrow(pladb_pdiff_introns_tab)

pladb_fstat_introns_tab <- pladb_fstat_introns[, c("GENE","EVENT","COORD","Fstat","deltapsi")] %>%
  filter(!is.na(EVENT)) %>%  # Remove rows where EVENT is NA
  arrange(desc(Fstat), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(pladb_fstat_introns_tab) <- 1:nrow(pladb_fstat_introns_tab)

pladb_fdr_introns_tab <- pladb_fdr_introns[, c("GENE","EVENT","COORD","FDR","deltapsi")] %>%
  filter(!is.na(EVENT), FDR<=0.05) %>%  # Remove rows where EVENT is NA
  arrange((FDR), desc(abs(deltapsi)))
rownames(pladb_fdr_introns_tab) <- 1:nrow(pladb_fdr_introns_tab)

# Find common EVENTs across all tables
common_events <- Reduce(intersect, list(pladb_pdiff_introns_tab$EVENT, pladb_fstat_introns_tab$EVENT, pladb_fdr_introns_tab$EVENT))

# Filter each table to include only rows with common EVENTs
filtered_introns_diff_tab <- pladb_pdiff_introns_tab %>% filter(EVENT %in% common_events)
filtered_introns_fstat_tab <- pladb_fstat_introns_tab %>% filter(EVENT %in% common_events)
filtered_introns_fdr_tab <- pladb_fdr_introns_tab %>% filter(EVENT %in% common_events)


# Select only unique columns from each table
filtered_introns_diff_tab <- filtered_introns_diff_tab %>% select(GENE,EVENT,COORD, deltapsi,Pdiff) # Extra column is Pdiff
filtered_introns_fstat_tab <- filtered_introns_fstat_tab %>% select(EVENT, Fstat) # Extra column is Fstat
filtered_introns_fdr_tab <- filtered_introns_fdr_tab %>% select(EVENT, FDR) # Extra column is FDR

# Merge the data frames by EVENT
pladb_intron_combination_tab <- filtered_introns_diff_tab %>%
  inner_join(filtered_introns_fstat_tab, by = "EVENT") %>%
  inner_join(filtered_introns_fdr_tab, by = "EVENT") %>%
  distinct() %>%
  arrange(desc(abs(deltapsi)), desc(Fstat),desc(FDR), desc(Pdiff))

# Render DataTable with enhancements
datatable(
  pladb_intron_combination_tab,
  options = list(
    pageLength = 10,                      # Rows per page
    autoWidth = TRUE,                     # Adjust column widths automatically
    dom = 'Bfrtip',                       # Add buttons for export
    buttons = c("copy", "csv", "excel", "pdf"),  # Simplified button definitions
    columnDefs = list(
      list(targets = "_all", className = "dt-center"),    # Center-align all columns
      list(
        targets = 3,  # Highlight significant Pdiff
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && parseFloat(data) > 0.99) {
               return '<span style=\"color: green; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      ),
      list(
        targets = 4,  # Highlight large deltapsi values
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && Math.abs(parseFloat(data)) > 0.5) {
               return '<span style=\"color: red; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      )
    )
  ),
  rownames = FALSE,                       # Disable row names
  extensions = "Buttons"                  # Enable export buttons
) %>%
  formatStyle(
    columns = "Pdiff",
    color = styleInterval(0.99, c("black", "green")),
    fontWeight = styleEqual(0.99, "bold")  # Ensure consistent styling
  ) %>%
  formatStyle(
  columns = "deltapsi",
  color = styleInterval(c(-0.5, 0.5), c("black", "red", "black")),  # Three colors for two intervals
  fontWeight = styleInterval(c(-0.5, 0.5), c("normal", "bold", "normal"))  # Three weights for two intervals
)

write.csv(pladb_intron_combination_tab,
          file = paste0("pladb_intron_combination_table_", format(Sys.Date(), "%Y%m%d"), ".csv"),
          row.names = FALSE)

3.2.3.2 PladB 1mM vs Control

pladb_pdiff_introns_tab <- pladb_1mm_control_pdiff_introns[, c("GENE", "EVENT", "COORD", "Pdiff", "deltapsi")] %>%
  filter(!is.na(EVENT), Pdiff >= 0.95) %>%  # Remove rows where EVENT is NA and Pdiff >= 0.95
  arrange(desc(Pdiff), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(pladb_pdiff_introns_tab) <- 1:nrow(pladb_pdiff_introns_tab)

pladb_fstat_introns_tab <- pladb_1mm_control_fstat_introns[, c("GENE","EVENT","COORD","Fstat","deltapsi")] %>%
  filter(!is.na(EVENT)) %>%  # Remove rows where EVENT is NA
  arrange(desc(Fstat), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(pladb_fstat_introns_tab) <- 1:nrow(pladb_fstat_introns_tab)

pladb_fdr_introns_tab <- pladb_1mm_control_fdr_introns[, c("GENE","EVENT","COORD","FDR","deltapsi")] %>%
  filter(!is.na(EVENT), FDR<=0.05) %>%  # Remove rows where EVENT is NA
  arrange((FDR), desc(abs(deltapsi)))
rownames(pladb_fdr_introns_tab) <- 1:nrow(pladb_fdr_introns_tab)

# Find common EVENTs across all tables
common_events <- Reduce(intersect, list(pladb_pdiff_introns_tab$EVENT, pladb_fstat_introns_tab$EVENT, pladb_fdr_introns_tab$EVENT))

# Filter each table to include only rows with common EVENTs
filtered_introns_diff_tab <- pladb_pdiff_introns_tab %>% filter(EVENT %in% common_events)
filtered_introns_fstat_tab <- pladb_fstat_introns_tab %>% filter(EVENT %in% common_events)
filtered_introns_fdr_tab <- pladb_fdr_introns_tab %>% filter(EVENT %in% common_events)


# Select only unique columns from each table
filtered_introns_diff_tab <- filtered_introns_diff_tab %>% select(GENE,EVENT,COORD, deltapsi,Pdiff) # Extra column is Pdiff
filtered_introns_fstat_tab <- filtered_introns_fstat_tab %>% select(EVENT, Fstat) # Extra column is Fstat
filtered_introns_fdr_tab <- filtered_introns_fdr_tab %>% select(EVENT, FDR) # Extra column is FDR

# Merge the data frames by EVENT
pladb_1mm_control_intron_combination_tab <- filtered_introns_diff_tab %>%
  inner_join(filtered_introns_fstat_tab, by = "EVENT") %>%
  inner_join(filtered_introns_fdr_tab, by = "EVENT") %>%
  distinct() %>%
  arrange(desc(abs(deltapsi)), desc(Fstat),desc(FDR), desc(Pdiff))

# Find if events appear in the control to 10mM pairwise comparison
pladb_1mm_control_intron_combination_tab$appears_in_10mm_control <- pladb_1mm_control_intron_combination_tab$EVENT %in% pladb_intron_combination_tab$EVENT

# Render DataTable with enhancements
datatable(
  pladb_1mm_control_intron_combination_tab,
  options = list(
    pageLength = 10,                      # Rows per page
    autoWidth = TRUE,                     # Adjust column widths automatically
    dom = 'Bfrtip',                       # Add buttons for export
    buttons = c("copy", "csv", "excel", "pdf"),  # Simplified button definitions
    columnDefs = list(
      list(targets = "_all", className = "dt-center"),    # Center-align all columns
      list(
        targets = 3,  # Highlight significant Pdiff
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && parseFloat(data) > 0.99) {
               return '<span style=\"color: green; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      ),
      list(
        targets = 4,  # Highlight large deltapsi values
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && Math.abs(parseFloat(data)) > 0.5) {
               return '<span style=\"color: red; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      )
    )
  ),
  rownames = FALSE,                       # Disable row names
  extensions = "Buttons"                  # Enable export buttons
) %>%
  formatStyle(
    columns = "Pdiff",
    color = styleInterval(0.99, c("black", "green")),
    fontWeight = styleEqual(0.99, "bold")  # Ensure consistent styling
  ) %>%
  formatStyle(
  columns = "deltapsi",
  color = styleInterval(c(-0.5, 0.5), c("black", "red", "black")),  # Three colors for two intervals
  fontWeight = styleInterval(c(-0.5, 0.5), c("normal", "bold", "normal"))  # Three weights for two intervals
)

write.csv(pladb_1mm_control_intron_combination_tab,
          file = paste0("pladb_1mm_control_intron_combination_table_", format(Sys.Date(), "%Y%m%d"), ".csv"),
          row.names = FALSE)

3.2.3.3 PladB 10mM vs 1mM

pladb_pdiff_introns_tab <- pladb_1mm_10mm_pdiff_introns[, c("GENE", "EVENT", "COORD", "Pdiff", "deltapsi")] %>%
  filter(!is.na(EVENT), Pdiff >= 0.95) %>%  # Remove rows where EVENT is NA and Pdiff >= 0.95
  arrange(desc(Pdiff), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(pladb_pdiff_introns_tab) <- 1:nrow(pladb_pdiff_introns_tab)

pladb_fstat_introns_tab <- pladb_1mm_10mm_fstat_introns[, c("GENE","EVENT","COORD","Fstat","deltapsi")] %>%
  filter(!is.na(EVENT)) %>%  # Remove rows where EVENT is NA
  arrange(desc(Fstat), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(pladb_fstat_introns_tab) <- 1:nrow(pladb_fstat_introns_tab)

pladb_fdr_introns_tab <- pladb_1mm_10mm_fdr_introns[, c("GENE","EVENT","COORD","FDR","deltapsi")] %>%
  filter(!is.na(EVENT), FDR<=0.05) %>%  # Remove rows where EVENT is NA
  arrange((FDR), desc(abs(deltapsi)))
rownames(pladb_fdr_introns_tab) <- 1:nrow(pladb_fdr_introns_tab)

# Find common EVENTs across all tables
common_events <- Reduce(intersect, list(pladb_pdiff_introns_tab$EVENT, pladb_fstat_introns_tab$EVENT, pladb_fdr_introns_tab$EVENT))

# Filter each table to include only rows with common EVENTs
filtered_introns_diff_tab <- pladb_pdiff_introns_tab %>% filter(EVENT %in% common_events)
filtered_introns_fstat_tab <- pladb_fstat_introns_tab %>% filter(EVENT %in% common_events)
filtered_introns_fdr_tab <- pladb_fdr_introns_tab %>% filter(EVENT %in% common_events)


# Select only unique columns from each table
filtered_introns_diff_tab <- filtered_introns_diff_tab %>% select(GENE,EVENT,COORD, deltapsi,Pdiff) # Extra column is Pdiff
filtered_introns_fstat_tab <- filtered_introns_fstat_tab %>% select(EVENT, Fstat) # Extra column is Fstat
filtered_introns_fdr_tab <- filtered_introns_fdr_tab %>% select(EVENT, FDR) # Extra column is FDR

# Merge the data frames by EVENT
pladb_10mm_1mm_intron_combination_tab <- filtered_introns_diff_tab %>%
  inner_join(filtered_introns_fstat_tab, by = "EVENT") %>%
  inner_join(filtered_introns_fdr_tab, by = "EVENT") %>%
  distinct() %>%
  arrange(desc(abs(deltapsi)), desc(Fstat),desc(FDR), desc(Pdiff))

# Find if events appear in the control to 10mM pairwise comparison
pladb_10mm_1mm_intron_combination_tab$appears_in_10mm_control <- pladb_10mm_1mm_intron_combination_tab$EVENT %in% pladb_intron_combination_tab$EVENT

# Render DataTable with enhancements
datatable(
  pladb_10mm_1mm_intron_combination_tab,
  options = list(
    pageLength = 10,                      # Rows per page
    autoWidth = TRUE,                     # Adjust column widths automatically
    dom = 'Bfrtip',                       # Add buttons for export
    buttons = c("copy", "csv", "excel", "pdf"),  # Simplified button definitions
    columnDefs = list(
      list(targets = "_all", className = "dt-center"),    # Center-align all columns
      list(
        targets = 3,  # Highlight significant Pdiff
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && parseFloat(data) > 0.99) {
               return '<span style=\"color: green; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      ),
      list(
        targets = 4,  # Highlight large deltapsi values
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && Math.abs(parseFloat(data)) > 0.5) {
               return '<span style=\"color: red; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      )
    )
  ),
  rownames = FALSE,                       # Disable row names
  extensions = "Buttons"                  # Enable export buttons
) %>%
  formatStyle(
    columns = "Pdiff",
    color = styleInterval(0.99, c("black", "green")),
    fontWeight = styleEqual(0.99, "bold")  # Ensure consistent styling
  ) %>%
  formatStyle(
  columns = "deltapsi",
  color = styleInterval(c(-0.5, 0.5), c("black", "red", "black")),  # Three colors for two intervals
  fontWeight = styleInterval(c(-0.5, 0.5), c("normal", "bold", "normal"))  # Three weights for two intervals
)

write.csv(pladb_10mm_1mm_intron_combination_tab,
          file = paste0("pladb_10mm_1mm_intron_combination_table_", format(Sys.Date(), "%Y%m%d"), ".csv"),
          row.names = FALSE)

3.2.3.4 FMN2DKO

fmndko_pdiff_introns_tab <- fmndko_pdiff_introns[, c("GENE", "EVENT", "COORD", "Pdiff", "deltapsi")] %>%
  filter(!is.na(EVENT), Pdiff >= 0.95) %>%  # Remove rows where EVENT is NA and Pdiff >= 0.95
  arrange(desc(Pdiff), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(fmndko_pdiff_introns_tab) <- 1:nrow(fmndko_pdiff_introns_tab)

fmndko_fstat_introns_tab <- fmndko_fstat_introns[, c("GENE","EVENT","COORD","Fstat","deltapsi")] %>%
  filter(!is.na(EVENT)) %>%  # Remove rows where EVENT is NA
  arrange(desc(Fstat), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(fmndko_fstat_introns_tab) <- 1:nrow(fmndko_fstat_introns_tab)

fmndko_fdr_introns_tab <- fmndko_fdr_introns[, c("GENE","EVENT","COORD","FDR","deltapsi")] %>%
  filter(!is.na(EVENT), FDR<=0.05) %>%  # Remove rows where EVENT is NA
  arrange((FDR), desc(abs(deltapsi)))
rownames(fmndko_fdr_introns_tab) <- 1:nrow(fmndko_fdr_introns_tab)

# Find common EVENTs across all tables
common_events <- Reduce(intersect, list(fmndko_pdiff_introns_tab$EVENT, fmndko_fstat_introns_tab$EVENT, fmndko_fdr_introns_tab$EVENT))

# Filter each table to include only rows with common EVENTs
filtered_intron_diff_tab <- fmndko_pdiff_introns_tab %>% filter(EVENT %in% common_events)
filtered_intron_fstat_tab <- fmndko_fstat_introns_tab %>% filter(EVENT %in% common_events)
filtered_intron_fdr_tab <- fmndko_fdr_introns_tab %>% filter(EVENT %in% common_events)


# Select only unique columns from each table
filtered_intron_diff_tab <- filtered_intron_diff_tab %>% select(GENE,EVENT,COORD, deltapsi,Pdiff) # Extra column is Pdiff
filtered_intron_fstat_tab <- filtered_intron_fstat_tab %>% select(EVENT, Fstat) # Extra column is Fstat
filtered_intron_fdr_tab <- filtered_intron_fdr_tab %>% select(EVENT, FDR) # Extra column is FDR

# Merge the data frames by EVENT
fmndko_intron_combination_tab <- filtered_intron_diff_tab %>%
  inner_join(filtered_intron_fstat_tab, by = "EVENT") %>%
  inner_join(filtered_intron_fdr_tab, by = "EVENT") %>%
  distinct() %>%
  arrange(desc(abs(deltapsi)), desc(Fstat),desc(FDR), desc(Pdiff))

# Render DataTable with enhancements
datatable(
  fmndko_intron_combination_tab,
  options = list(
    pageLength = 10,                      # Rows per page
    autoWidth = TRUE,                     # Adjust column widths automatically
    dom = 'Bfrtip',                       # Add buttons for export
    buttons = c("copy", "csv", "excel", "pdf"),  # Simplified button definitions
    columnDefs = list(
      list(targets = "_all", className = "dt-center"),    # Center-align all columns
      list(
        targets = 3,  # Highlight significant Pdiff
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && parseFloat(data) > 0.99) {
               return '<span style=\"color: green; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      ),
      list(
        targets = 4,  # Highlight large deltapsi values
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && Math.abs(parseFloat(data)) > 0.5) {
               return '<span style=\"color: red; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      )
    )
  ),
  rownames = FALSE,                       # Disable row names
  extensions = "Buttons"                  # Enable export buttons
) %>%
  formatStyle(
    columns = "Pdiff",
    color = styleInterval(0.99, c("black", "green")),
    fontWeight = styleEqual(0.99, "bold")  # Ensure consistent styling
  ) %>%
  formatStyle(
  columns = "deltapsi",
  color = styleInterval(c(-0.5, 0.5), c("black", "red", "black")),  # Three colors for two intervals
  fontWeight = styleInterval(c(-0.5, 0.5), c("normal", "bold", "normal"))  # Three weights for two intervals
)

write.csv(fmndko_intron_combination_tab,
          file = paste0("fmndko_intron_combination_table_", format(Sys.Date(), "%Y%m%d"), ".csv"),
          row.names = FALSE)

3.2.3.5 Spire2

spire_pdiff_introns_tab <- spire_pdiff_introns[, c("GENE", "EVENT", "COORD", "Pdiff", "deltapsi")] %>%
  filter(!is.na(EVENT), Pdiff >= 0.95) %>%  # Remove rows where EVENT is NA and Pdiff >= 0.95
  arrange(desc(Pdiff), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(spire_pdiff_introns_tab) <- 1:nrow(spire_pdiff_introns_tab)

spire_fstat_introns_tab <- spire_fstat_introns[, c("GENE","EVENT","COORD","Fstat","deltapsi")] %>%
  filter(!is.na(EVENT)) %>%  # Remove rows where EVENT is NA
  arrange(desc(Fstat), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(spire_fstat_introns_tab) <- 1:nrow(spire_fstat_introns_tab)

spire_fdr_introns_tab <- spire_fdr_introns[, c("GENE","EVENT","COORD","FDR","deltapsi")] %>%
  filter(!is.na(EVENT), FDR<=0.05) %>%  # Remove rows where EVENT is NA
  arrange((FDR), desc(abs(deltapsi)))
rownames(spire_fdr_introns_tab) <- 1:nrow(spire_fdr_introns_tab)

# Find common EVENTs across all tables
common_events <- Reduce(intersect, list(spire_pdiff_introns_tab$EVENT, spire_fstat_introns_tab$EVENT, spire_fdr_introns_tab$EVENT))

# Filter each table to include only rows with common EVENTs
filtered_intron_diff_tab <- spire_pdiff_introns_tab %>% filter(EVENT %in% common_events)
filtered_intron_fstat_tab <- spire_fstat_introns_tab %>% filter(EVENT %in% common_events)
filtered_intron_fdr_tab <- spire_fdr_introns_tab %>% filter(EVENT %in% common_events)


# Select only unique columns from each table
filtered_intron_diff_tab <- filtered_intron_diff_tab %>% select(GENE,EVENT,COORD, deltapsi,Pdiff) # Extra column is Pdiff
filtered_intron_fstat_tab <- filtered_intron_fstat_tab %>% select(EVENT, Fstat) # Extra column is Fstat
filtered_intron_fdr_tab <- filtered_intron_fdr_tab %>% select(EVENT, FDR) # Extra column is FDR

# Merge the data frames by EVENT
spire_intron_combination_tab <- filtered_intron_diff_tab %>%
  inner_join(filtered_intron_fstat_tab, by = "EVENT") %>%
  inner_join(filtered_intron_fdr_tab, by = "EVENT") %>%
  distinct() %>%
  arrange(desc(abs(deltapsi)), desc(Fstat),desc(FDR), desc(Pdiff))

# Render DataTable with enhancements
datatable(
  spire_intron_combination_tab,
  options = list(
    pageLength = 10,                      # Rows per page
    autoWidth = TRUE,                     # Adjust column widths automatically
    dom = 'Bfrtip',                       # Add buttons for export
    buttons = c("copy", "csv", "excel", "pdf"),  # Simplified button definitions
    columnDefs = list(
      list(targets = "_all", className = "dt-center"),    # Center-align all columns
      list(
        targets = 3,  # Highlight significant Pdiff
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && parseFloat(data) > 0.99) {
               return '<span style=\"color: green; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      ),
      list(
        targets = 4,  # Highlight large deltapsi values
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && Math.abs(parseFloat(data)) > 0.5) {
               return '<span style=\"color: red; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      )
    )
  ),
  rownames = FALSE,                       # Disable row names
  extensions = "Buttons"                  # Enable export buttons
) %>%
  formatStyle(
    columns = "Pdiff",
    color = styleInterval(0.99, c("black", "green")),
    fontWeight = styleEqual(0.99, "bold")  # Ensure consistent styling
  ) %>%
  formatStyle(
  columns = "deltapsi",
  color = styleInterval(c(-0.5, 0.5), c("black", "red", "black")),  # Three colors for two intervals
  fontWeight = styleInterval(c(-0.5, 0.5), c("normal", "bold", "normal"))  # Three weights for two intervals
)

write.csv(spire_intron_combination_tab,
          file = paste0("spire_intron_combination_table_", format(Sys.Date(), "%Y%m%d"), ".csv"),
          row.names = FALSE)

3.3 Alternative Splicing List

3.3.1 Calculations

pladb_pdiff_alt <- prepareTableVolcano(
  psitable = pladb_alt$PSI,
  qualtable = pladb_alt$Qual,
  npoints = 500,
  colsA = colsGroupA_pladb,
  colsB = colsGroupB_pladb,
  labA = groupA_pladb,
  labB = groupB_pladb,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed=TRUE,
  CoverageWeight = FALSE)

pladb_1mm_control_pdiff_alt <- prepareTableVolcano(
  psitable = pladb_alt$PSI,
  qualtable = pladb_alt$Qual,
  npoints = 500,
  colsA = colsGroupA_pladb,
  colsB = colsGroupC_pladb,
  labA = groupA_pladb,
  labB = groupC_pladb,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed=TRUE,
  CoverageWeight = FALSE)

pladb_1mm_10mm_pdiff_alt <- prepareTableVolcano(
  psitable = pladb_alt$PSI,
  qualtable = pladb_alt$Qual,
  npoints = 500,
  colsA = colsGroupC_pladb,
  colsB = colsGroupB_pladb,
  labA = groupC_pladb,
  labB = groupB_pladb,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed=TRUE,
  CoverageWeight = FALSE)

# Prepare table for Tao
fmndko_pdiff_alt <- prepareTableVolcano(
  psitable = fmndko_alt$PSI,
  qualtable = fmndko_alt$Qual,
  npoints = 500,
  colsA = colsGroupA_fmndko,
  colsB = colsGroupB_fmndko,
  labA = groupA_fmndko,
  labB = groupB_fmndko,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed = TRUE,
  CoverageWeight = FALSE
)

# Prepare table for spire
spire_pdiff_alt <- prepareTableVolcano(
  psitable = spire_alt$PSI,
  qualtable = spire_alt$Qual,
  npoints = 500,
  colsA = colsGroupA_spire,
  colsB = colsGroupB_spire,
  labA = groupA_spire,
  labB = groupB_spire,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed = TRUE,
  CoverageWeight = FALSE
)

pladb_fstat_alt <- prepareTableVolcanoFstat(
  psitable = pladb_alt$PSI,
  qualtable = pladb_alt$Qual,
  npoints = 500,
  colsA = colsGroupA_pladb,
  colsB = colsGroupB_pladb,
  labA = groupA_pladb,
  labB = groupB_pladb,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed=TRUE,
  CoverageWeight = FALSE)

pladb_1mm_control_fstat_alt <- prepareTableVolcanoFstat(
  psitable = pladb_alt$PSI,
  qualtable = pladb_alt$Qual,
  npoints = 500,
  colsA = colsGroupA_pladb,
  colsB = colsGroupC_pladb,
  labA = groupA_pladb,
  labB = groupC_pladb,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed=TRUE,
  CoverageWeight = FALSE)

pladb_1mm_10mm_fstat_alt <- prepareTableVolcanoFstat(
  psitable = pladb_alt$PSI,
  qualtable = pladb_alt$Qual,
  npoints = 500,
  colsA = colsGroupC_pladb,
  colsB = colsGroupB_pladb,
  labA = groupC_pladb,
  labB = groupB_pladb,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed=TRUE,
  CoverageWeight = FALSE)

# Prepare table for fmndko
fmndko_fstat_alt <- prepareTableVolcanoFstat(
  psitable = fmndko_alt$PSI,
  qualtable = fmndko_alt$Qual,
  npoints = 500,
  colsA = colsGroupA_fmndko,
  colsB = colsGroupB_fmndko,
  labA = groupA_fmndko,
  labB = groupB_fmndko,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed = TRUE,
  CoverageWeight = FALSE
)


# Prepare table for spire
spire_fstat_alt <- prepareTableVolcanoFstat(
  psitable = spire_alt$PSI,
  qualtable = spire_alt$Qual,
  npoints = 500,
  colsA = colsGroupA_spire,
  colsB = colsGroupB_spire,
  labA = groupA_spire,
  labB = groupB_spire,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  seed = TRUE,
  CoverageWeight = FALSE
)

pladb_fdr_alt <- prepareTableVolcanoFDR(
  psitable = pladb_alt$PSI,
  qualtable = pladb_alt$Qual,
  npoints = 500,
  colsA = colsGroupA_pladb,
  colsB = colsGroupB_pladb,
  labA = groupA_pladb,
  labB = groupB_pladb,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  nsim = 100,
  seed=TRUE,
  CoverageWeight = FALSE)

pladb_1mm_control_fdr_alt <- prepareTableVolcanoFDR(
  psitable = pladb_alt$PSI,
  qualtable = pladb_alt$Qual,
  npoints = 500,
  colsA = colsGroupA_pladb,
  colsB = colsGroupC_pladb,
  labA = groupA_pladb,
  labB = groupC_pladb,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  nsim = 100,
  seed=TRUE,
  CoverageWeight = FALSE)


pladb_1mm_10mm_fdr_alt <- prepareTableVolcanoFDR(
  psitable = pladb_alt$PSI,
  qualtable = pladb_alt$Qual,
  npoints = 500,
  colsA = colsGroupC_pladb,
  colsB = colsGroupB_pladb,
  labA = groupC_pladb,
  labB = groupB_pladb,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  nsim = 100,
  seed=TRUE,
  CoverageWeight = FALSE)


# Prepare table for fmndko
fmndko_fdr_alt <- prepareTableVolcanoFDR(
  psitable = fmndko_alt$PSI,
  qualtable = fmndko_alt$Qual,
  npoints = 500,
  colsA = colsGroupA_fmndko,
  colsB = colsGroupB_fmndko,
  labA = groupA_fmndko,
  labB = groupB_fmndko,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  nsim = 100,
  seed = TRUE,
  CoverageWeight = FALSE
)

# Prepare table for spire
spire_fdr_alt <- prepareTableVolcanoFDR(
  psitable = spire_alt$PSI,
  qualtable = spire_alt$Qual,
  npoints = 500,
  colsA = colsGroupA_spire,
  colsB = colsGroupB_spire,
  labA = groupA_spire,
  labB = groupB_spire,
  basalColor = "#89C0AE",
  interestColor = "#E69A9C",
  maxDevTable = maxDevSimulationN100,
  nsim = 100,
  seed = TRUE,
  CoverageWeight = FALSE
)

3.3.2 FDR Volcano Plots

3.3.2.1 PladB 10mM vs Control

plotVolcanoFDR(betasTable =filter(pladb_fdr_alt,!is.na(EVENT)),
                            labA = groupA_pladb,
                            labB = groupB_pladb,
                            basalColor = "#89C0AE",
                            interestColor = "#E69A9C") +

theme(
  plot.title = element_text(size = 10),  # Title font size
  axis.title = element_text(size = 10),  # Axis title font size
  axis.text = element_text(size = 10),   # Axis text font size
  legend.text = element_text(size = 10), # Legend text font size
  legend.title = element_text(size = 10) # Legend title font size
)

3.3.2.2 PladB 1mM vs Control

plotVolcanoFDR(betasTable =filter(pladb_1mm_control_fdr_alt,!is.na(EVENT)),
                            labA = groupA_pladb,
                            labB = groupC_pladb,
                            basalColor = "#89C0AE",
                            interestColor = "#E69A9C") +

theme(
  plot.title = element_text(size = 10),  # Title font size
  axis.title = element_text(size = 10),  # Axis title font size
  axis.text = element_text(size = 10),   # Axis text font size
  legend.text = element_text(size = 10), # Legend text font size
  legend.title = element_text(size = 10) # Legend title font size
)

3.3.2.3 PladB 10mM vs 1mM

plotVolcanoFDR(betasTable =filter(pladb_1mm_10mm_fdr_alt,!is.na(EVENT)),
                            labA = groupC_pladb,
                            labB = groupB_pladb,
                            basalColor = "#89C0AE",
                            interestColor = "#E69A9C") +

theme(
  plot.title = element_text(size = 10),  # Title font size
  axis.title = element_text(size = 10),  # Axis title font size
  axis.text = element_text(size = 10),   # Axis text font size
  legend.text = element_text(size = 10), # Legend text font size
  legend.title = element_text(size = 10) # Legend title font size
)

3.3.2.4 FMN2DKO

plotVolcanoFDR(betasTable =filter(fmndko_fdr_alt,!is.na(EVENT)),
                            labA = groupA_fmndko,
                            labB = groupB_fmndko,
                            basalColor = "#89C0AE",
                            interestColor = "#E69A9C") +
theme(
  plot.title = element_text(size = 10),  # Title font size
  axis.title = element_text(size = 10),  # Axis title font size
  axis.text = element_text(size = 10),   # Axis text font size
  legend.text = element_text(size = 10), # Legend text font size
  legend.title = element_text(size = 10) # Legend title font size
)

3.3.2.5 Spire2

plotVolcanoFDR(betasTable =filter(spire_fdr_alt,!is.na(EVENT)),
                            labA = groupA_spire,
                            labB = groupB_spire,
                            basalColor = "#89C0AE",
                            interestColor = "#E69A9C") +
theme(
  plot.title = element_text(size = 10),  # Title font size
  axis.title = element_text(size = 10),  # Axis title font size
  axis.text = element_text(size = 10),   # Axis text font size
  legend.text = element_text(size = 10), # Legend text font size
  legend.title = element_text(size = 10) # Legend title font size
)

3.3.3 Individual Combination Tables

3.3.3.1 PladB 10mM vs control

pladb_pdiff_alt_tab <- pladb_pdiff_alt[, c("GENE", "EVENT", "COORD", "Pdiff", "deltapsi")] %>%
  filter(!is.na(EVENT), Pdiff >= 0.95) %>%  # Remove rows where EVENT is NA and Pdiff >= 0.95
  arrange(desc(Pdiff), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(pladb_pdiff_alt_tab) <- 1:nrow(pladb_pdiff_alt_tab)

pladb_fstat_alt_tab <- pladb_fstat_alt[, c("GENE","EVENT","COORD","Fstat","deltapsi")] %>%
  filter(!is.na(EVENT)) %>%  # Remove rows where EVENT is NA
  arrange(desc(Fstat), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(pladb_fstat_alt_tab) <- 1:nrow(pladb_fstat_alt_tab)

pladb_fdr_alt_tab <- pladb_fdr_alt[, c("GENE","EVENT","COORD","FDR","deltapsi")] %>%
  filter(!is.na(EVENT), FDR<=0.05) %>%  # Remove rows where EVENT is NA
  arrange((FDR), desc(abs(deltapsi)))
rownames(pladb_fdr_alt_tab) <- 1:nrow(pladb_fdr_alt_tab)

# Find common EVENTs across all tables
common_events <- Reduce(intersect, list(pladb_pdiff_alt_tab$EVENT, pladb_fstat_alt_tab$EVENT, pladb_fdr_alt_tab$EVENT))

# Filter each table to include only rows with common EVENTs
filtered_alt_diff_tab <- pladb_pdiff_alt_tab %>% filter(EVENT %in% common_events)
filtered_alt_fstat_tab <- pladb_fstat_alt_tab %>% filter(EVENT %in% common_events)
filtered_alt_fdr_tab <- pladb_fdr_alt_tab %>% filter(EVENT %in% common_events)


# Select only unique columns from each table
filtered_alt_diff_tab <- filtered_alt_diff_tab %>% select(GENE,EVENT,COORD, deltapsi,Pdiff) # Extra column is Pdiff
filtered_alt_fstat_tab <- filtered_alt_fstat_tab %>% select(EVENT, Fstat) # Extra column is Fstat
filtered_alt_fdr_tab <- filtered_alt_fdr_tab %>% select(EVENT, FDR) # Extra column is FDR

# Merge the data frames by EVENT
pladb_alt_combination_tab <- filtered_alt_diff_tab %>%
  inner_join(filtered_alt_fstat_tab, by = "EVENT") %>%
  inner_join(filtered_alt_fdr_tab, by = "EVENT") %>%
  distinct() %>%
  arrange(desc(abs(deltapsi)), desc(Fstat),desc(FDR), desc(Pdiff))

# Render DataTable with enhancements
datatable(
  pladb_alt_combination_tab,
  options = list(
    pageLength = 10,                      # Rows per page
    autoWidth = TRUE,                     # Adjust column widths automatically
    dom = 'Bfrtip',                       # Add buttons for export
    buttons = c("copy", "csv", "excel", "pdf"),  # Simplified button definitions
    columnDefs = list(
      list(targets = "_all", className = "dt-center"),    # Center-align all columns
      list(
        targets = 3,  # Highlight significant Pdiff
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && parseFloat(data) > 0.99) {
               return '<span style=\"color: green; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      ),
      list(
        targets = 4,  # Highlight large deltapsi values
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && Math.abs(parseFloat(data)) > 0.5) {
               return '<span style=\"color: red; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      )
    )
  ),
  rownames = FALSE,                       # Disable row names
  extensions = "Buttons"                  # Enable export buttons
) %>%
  formatStyle(
    columns = "Pdiff",
    color = styleInterval(0.99, c("black", "green")),
    fontWeight = styleEqual(0.99, "bold")  # Ensure consistent styling
  ) %>%
  formatStyle(
  columns = "deltapsi",
  color = styleInterval(c(-0.5, 0.5), c("black", "red", "black")),  # Three colors for two intervals
  fontWeight = styleInterval(c(-0.5, 0.5), c("normal", "bold", "normal"))  # Three weights for two intervals
)

write.csv(pladb_alt_combination_tab,
          file = paste0("pladb_alt_combination_table_", format(Sys.Date(), "%Y%m%d"), ".csv"),
          row.names = FALSE)

3.3.3.2 PladB 1mM vs Control

pladb_pdiff_alt_tab <- pladb_1mm_control_pdiff_alt[, c("GENE", "EVENT", "COORD", "Pdiff", "deltapsi")] %>%
  filter(!is.na(EVENT), Pdiff >= 0.95) %>%  # Remove rows where EVENT is NA and Pdiff >= 0.95
  arrange(desc(Pdiff), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(pladb_pdiff_alt_tab) <- 1:nrow(pladb_pdiff_alt_tab)

pladb_fstat_alt_tab <- pladb_1mm_control_fstat_alt[, c("GENE","EVENT","COORD","Fstat","deltapsi")] %>%
  filter(!is.na(EVENT)) %>%  # Remove rows where EVENT is NA
  arrange(desc(Fstat), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(pladb_fstat_alt_tab) <- 1:nrow(pladb_fstat_alt_tab)

pladb_fdr_alt_tab <- pladb_1mm_control_fdr_alt[, c("GENE","EVENT","COORD","FDR","deltapsi")] %>%
  filter(!is.na(EVENT), FDR<=0.05) %>%  # Remove rows where EVENT is NA
  arrange((FDR), desc(abs(deltapsi)))
rownames(pladb_fdr_alt_tab) <- 1:nrow(pladb_fdr_alt_tab)

# Find common EVENTs across all tables
common_events <- Reduce(intersect, list(pladb_pdiff_alt_tab$EVENT, pladb_fstat_alt_tab$EVENT, pladb_fdr_alt_tab$EVENT))

# Filter each table to include only rows with common EVENTs
filtered_alt_diff_tab <- pladb_pdiff_alt_tab %>% filter(EVENT %in% common_events)
filtered_alt_fstat_tab <- pladb_fstat_alt_tab %>% filter(EVENT %in% common_events)
filtered_alt_fdr_tab <- pladb_fdr_alt_tab %>% filter(EVENT %in% common_events)


# Select only unique columns from each table
filtered_alt_diff_tab <- filtered_alt_diff_tab %>% select(GENE,EVENT,COORD, deltapsi,Pdiff) # Extra column is Pdiff
filtered_alt_fstat_tab <- filtered_alt_fstat_tab %>% select(EVENT, Fstat) # Extra column is Fstat
filtered_alt_fdr_tab <- filtered_alt_fdr_tab %>% select(EVENT, FDR) # Extra column is FDR

# Merge the data frames by EVENT
pladb_1mm_control_alt_combination_tab <- filtered_alt_diff_tab %>%
  inner_join(filtered_alt_fstat_tab, by = "EVENT") %>%
  inner_join(filtered_alt_fdr_tab, by = "EVENT") %>%
  distinct() %>%
  arrange(desc(abs(deltapsi)), desc(Fstat),desc(FDR), desc(Pdiff))

# Find if events appear in the control to 10mM pairwise comparison
pladb_1mm_control_alt_combination_tab$appears_in_10mm_control <- pladb_1mm_control_alt_combination_tab$EVENT %in% pladb_alt_combination_tab$EVENT

# Render DataTable with enhancements
datatable(
  pladb_1mm_control_alt_combination_tab,
  options = list(
    pageLength = 10,                      # Rows per page
    autoWidth = TRUE,                     # Adjust column widths automatically
    dom = 'Bfrtip',                       # Add buttons for export
    buttons = c("copy", "csv", "excel", "pdf"),  # Simplified button definitions
    columnDefs = list(
      list(targets = "_all", className = "dt-center"),    # Center-align all columns
      list(
        targets = 3,  # Highlight significant Pdiff
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && parseFloat(data) > 0.99) {
               return '<span style=\"color: green; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      ),
      list(
        targets = 4,  # Highlight large deltapsi values
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && Math.abs(parseFloat(data)) > 0.5) {
               return '<span style=\"color: red; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      )
    )
  ),
  rownames = FALSE,                       # Disable row names
  extensions = "Buttons"                  # Enable export buttons
) %>%
  formatStyle(
    columns = "Pdiff",
    color = styleInterval(0.99, c("black", "green")),
    fontWeight = styleEqual(0.99, "bold")  # Ensure consistent styling
  ) %>%
  formatStyle(
  columns = "deltapsi",
  color = styleInterval(c(-0.5, 0.5), c("black", "red", "black")),  # Three colors for two intervals
  fontWeight = styleInterval(c(-0.5, 0.5), c("normal", "bold", "normal"))  # Three weights for two intervals
)

write.csv(pladb_1mm_control_alt_combination_tab,
          file = paste0("pladb_1mm_control_alt_combination_table_", format(Sys.Date(), "%Y%m%d"), ".csv"),
          row.names = FALSE)

3.3.3.3 PladB 10mM vs 1mM

pladb_pdiff_alt_tab <- pladb_1mm_10mm_pdiff_alt[, c("GENE", "EVENT", "COORD", "Pdiff", "deltapsi")] %>%
  filter(!is.na(EVENT), Pdiff >= 0.95) %>%  # Remove rows where EVENT is NA and Pdiff >= 0.95
  arrange(desc(Pdiff), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(pladb_pdiff_alt_tab) <- 1:nrow(pladb_pdiff_alt_tab)

pladb_fstat_alt_tab <- pladb_1mm_10mm_fstat_alt[, c("GENE","EVENT","COORD","Fstat","deltapsi")] %>%
  filter(!is.na(EVENT)) %>%  # Remove rows where EVENT is NA
  arrange(desc(Fstat), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(pladb_fstat_alt_tab) <- 1:nrow(pladb_fstat_alt_tab)

pladb_fdr_alt_tab <- pladb_1mm_10mm_fdr_alt[, c("GENE","EVENT","COORD","FDR","deltapsi")] %>%
  filter(!is.na(EVENT), FDR<=0.05) %>%  # Remove rows where EVENT is NA
  arrange((FDR), desc(abs(deltapsi)))
rownames(pladb_fdr_alt_tab) <- 1:nrow(pladb_fdr_alt_tab)

# Find common EVENTs across all tables
common_events <- Reduce(intersect, list(pladb_pdiff_alt_tab$EVENT, pladb_fstat_alt_tab$EVENT, pladb_fdr_alt_tab$EVENT))

# Filter each table to include only rows with common EVENTs
filtered_alt_diff_tab <- pladb_pdiff_alt_tab %>% filter(EVENT %in% common_events)
filtered_alt_fstat_tab <- pladb_fstat_alt_tab %>% filter(EVENT %in% common_events)
filtered_alt_fdr_tab <- pladb_fdr_alt_tab %>% filter(EVENT %in% common_events)


# Select only unique columns from each table
filtered_alt_diff_tab <- filtered_alt_diff_tab %>% select(GENE,EVENT,COORD, deltapsi,Pdiff) # Extra column is Pdiff
filtered_alt_fstat_tab <- filtered_alt_fstat_tab %>% select(EVENT, Fstat) # Extra column is Fstat
filtered_alt_fdr_tab <- filtered_alt_fdr_tab %>% select(EVENT, FDR) # Extra column is FDR

# Merge the data frames by EVENT
pladb_10mm_1mm_alt_combination_tab <- filtered_alt_diff_tab %>%
  inner_join(filtered_alt_fstat_tab, by = "EVENT") %>%
  inner_join(filtered_alt_fdr_tab, by = "EVENT") %>%
  distinct() %>%
  arrange(desc(abs(deltapsi)), desc(Fstat),desc(FDR), desc(Pdiff))

pladb_10mm_1mm_alt_combination_tab$appears_in_10mm_control <- pladb_10mm_1mm_alt_combination_tab$EVENT %in% pladb_alt_combination_tab$EVENT

# Render DataTable with enhancements
datatable(
  pladb_10mm_1mm_alt_combination_tab,
  options = list(
    pageLength = 10,                      # Rows per page
    autoWidth = TRUE,                     # Adjust column widths automatically
    dom = 'Bfrtip',                       # Add buttons for export
    buttons = c("copy", "csv", "excel", "pdf"),  # Simplified button definitions
    columnDefs = list(
      list(targets = "_all", className = "dt-center"),    # Center-align all columns
      list(
        targets = 3,  # Highlight significant Pdiff
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && parseFloat(data) > 0.99) {
               return '<span style=\"color: green; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      ),
      list(
        targets = 4,  # Highlight large deltapsi values
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && Math.abs(parseFloat(data)) > 0.5) {
               return '<span style=\"color: red; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      )
    )
  ),
  rownames = FALSE,                       # Disable row names
  extensions = "Buttons"                  # Enable export buttons
) %>%
  formatStyle(
    columns = "Pdiff",
    color = styleInterval(0.99, c("black", "green")),
    fontWeight = styleEqual(0.99, "bold")  # Ensure consistent styling
  ) %>%
  formatStyle(
  columns = "deltapsi",
  color = styleInterval(c(-0.5, 0.5), c("black", "red", "black")),  # Three colors for two intervals
  fontWeight = styleInterval(c(-0.5, 0.5), c("normal", "bold", "normal"))  # Three weights for two intervals
)

write.csv(pladb_10mm_1mm_alt_combination_tab,
          file = paste0("pladb_10mm_1mm_alt_combination_table_", format(Sys.Date(), "%Y%m%d"), ".csv"),
          row.names = FALSE)

3.3.3.4 FMN2DKO

fmndko_pdiff_alt_tab <- fmndko_pdiff_alt[, c("GENE", "EVENT", "COORD", "Pdiff", "deltapsi")] %>%
  filter(!is.na(EVENT), Pdiff >= 0.95) %>%  # Remove rows where EVENT is NA and Pdiff >= 0.95
  arrange(desc(Pdiff), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(fmndko_pdiff_alt_tab) <- 1:nrow(fmndko_pdiff_alt_tab)

fmndko_fstat_alt_tab <- fmndko_fstat_alt[, c("GENE","EVENT","COORD","Fstat","deltapsi")] %>%
  filter(!is.na(EVENT)) %>%  # Remove rows where EVENT is NA
  arrange(desc(Fstat), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(fmndko_fstat_alt_tab) <- 1:nrow(fmndko_fstat_alt_tab)

fmndko_fdr_alt_tab <- fmndko_fdr_alt[, c("GENE","EVENT","COORD","FDR","deltapsi")] %>%
  filter(!is.na(EVENT), FDR<=0.05) %>%  # Remove rows where EVENT is NA
  arrange((FDR), desc(abs(deltapsi)))
rownames(fmndko_fdr_alt_tab) <- 1:nrow(fmndko_fdr_alt_tab)

# Find common EVENTs across all tables
common_events <- Reduce(intersect, list(fmndko_pdiff_alt_tab$EVENT, fmndko_fstat_alt_tab$EVENT, fmndko_fdr_alt_tab$EVENT))

# Filter each table to include only rows with common EVENTs
filtered_alt_diff_tab <- fmndko_pdiff_alt_tab %>% filter(EVENT %in% common_events)
filtered_alt_fstat_tab <- fmndko_fstat_alt_tab %>% filter(EVENT %in% common_events)
filtered_alt_fdr_tab <- fmndko_fdr_alt_tab %>% filter(EVENT %in% common_events)


# Select only unique columns from each table
filtered_alt_diff_tab <- filtered_alt_diff_tab %>% select(GENE,EVENT,COORD, deltapsi,Pdiff) # Extra column is Pdiff
filtered_alt_fstat_tab <- filtered_alt_fstat_tab %>% select(EVENT, Fstat) # Extra column is Fstat
filtered_alt_fdr_tab <- filtered_alt_fdr_tab %>% select(EVENT, FDR) # Extra column is FDR

# Merge the data frames by EVENT
fmndko_alt_combination_tab <- filtered_alt_diff_tab %>%
  inner_join(filtered_alt_fstat_tab, by = "EVENT") %>%
  inner_join(filtered_alt_fdr_tab, by = "EVENT") %>%
  distinct() %>%
  arrange(desc(abs(deltapsi)), desc(Fstat),desc(FDR), desc(Pdiff))

# Render DataTable with enhancements
datatable(
  fmndko_alt_combination_tab,
  options = list(
    pageLength = 10,                      # Rows per page
    autoWidth = TRUE,                     # Adjust column widths automatically
    dom = 'Bfrtip',                       # Add buttons for export
    buttons = c("copy", "csv", "excel", "pdf"),  # Simplified button definitions
    columnDefs = list(
      list(targets = "_all", className = "dt-center"),    # Center-align all columns
      list(
        targets = 3,  # Highlight significant Pdiff
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && parseFloat(data) > 0.99) {
               return '<span style=\"color: green; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      ),
      list(
        targets = 4,  # Highlight large deltapsi values
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && Math.abs(parseFloat(data)) > 0.5) {
               return '<span style=\"color: red; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      )
    )
  ),
  rownames = FALSE,                       # Disable row names
  extensions = "Buttons"                  # Enable export buttons
) %>%
  formatStyle(
    columns = "Pdiff",
    color = styleInterval(0.99, c("black", "green")),
    fontWeight = styleEqual(0.99, "bold")  # Ensure consistent styling
  ) %>%
  formatStyle(
  columns = "deltapsi",
  color = styleInterval(c(-0.5, 0.5), c("black", "red", "black")),  # Three colors for two intervals
  fontWeight = styleInterval(c(-0.5, 0.5), c("normal", "bold", "normal"))  # Three weights for two intervals
)

write.csv(fmndko_alt_combination_tab,
          file = paste0("fmndko_alt_combination_table_", format(Sys.Date(), "%Y%m%d"), ".csv"),
          row.names = FALSE)

3.3.3.5 Spire2

spire_pdiff_alt_tab <- spire_pdiff_alt[, c("GENE", "EVENT", "COORD", "Pdiff", "deltapsi")] %>%
  filter(!is.na(EVENT), Pdiff >= 0.95) %>%  # Remove rows where EVENT is NA and Pdiff >= 0.95
  arrange(desc(Pdiff), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(spire_pdiff_alt_tab) <- 1:nrow(spire_pdiff_alt_tab)

spire_fstat_alt_tab <- spire_fstat_alt[, c("GENE","EVENT","COORD","Fstat","deltapsi")] %>%
  filter(!is.na(EVENT)) %>%  # Remove rows where EVENT is NA
  arrange(desc(Fstat), desc(abs(deltapsi)))  # Sort by Pdiff and abs(deltapsi)
rownames(spire_fstat_alt_tab) <- 1:nrow(spire_fstat_alt_tab)

spire_fdr_alt_tab <- spire_fdr_alt[, c("GENE","EVENT","COORD","FDR","deltapsi")] %>%
  filter(!is.na(EVENT), FDR<=0.05) %>%  # Remove rows where EVENT is NA
  arrange((FDR), desc(abs(deltapsi)))
rownames(spire_fdr_alt_tab) <- 1:nrow(spire_fdr_alt_tab)

# Find common EVENTs across all tables
common_events <- Reduce(intersect, list(spire_pdiff_alt_tab$EVENT, spire_fstat_alt_tab$EVENT, spire_fdr_alt_tab$EVENT))

# Filter each table to include only rows with common EVENTs
filtered_alt_diff_tab <- spire_pdiff_alt_tab %>% filter(EVENT %in% common_events)
filtered_alt_fstat_tab <- spire_fstat_alt_tab %>% filter(EVENT %in% common_events)
filtered_alt_fdr_tab <- spire_fdr_alt_tab %>% filter(EVENT %in% common_events)


# Select only unique columns from each table
filtered_alt_diff_tab <- filtered_alt_diff_tab %>% select(GENE,EVENT,COORD, deltapsi,Pdiff) # Extra column is Pdiff
filtered_alt_fstat_tab <- filtered_alt_fstat_tab %>% select(EVENT, Fstat) # Extra column is Fstat
filtered_alt_fdr_tab <- filtered_alt_fdr_tab %>% select(EVENT, FDR) # Extra column is FDR

# Merge the data frames by EVENT
spire_alt_combination_tab <- filtered_alt_diff_tab %>%
  inner_join(filtered_alt_fstat_tab, by = "EVENT") %>%
  inner_join(filtered_alt_fdr_tab, by = "EVENT") %>%
  distinct() %>%
  arrange(desc(abs(deltapsi)), desc(Fstat),desc(FDR), desc(Pdiff))

# Render DataTable with enhancements
datatable(
  spire_alt_combination_tab,
  options = list(
    pageLength = 10,                      # Rows per page
    autoWidth = TRUE,                     # Adjust column widths automatically
    dom = 'Bfrtip',                       # Add buttons for export
    buttons = c("copy", "csv", "excel", "pdf"),  # Simplified button definitions
    columnDefs = list(
      list(targets = "_all", className = "dt-center"),    # Center-align all columns
      list(
        targets = 3,  # Highlight significant Pdiff
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && parseFloat(data) > 0.99) {
               return '<span style=\"color: green; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      ),
      list(
        targets = 4,  # Highlight large deltapsi values
        render = JS(
          "function(data, type, row) {
             if (type === 'display' && Math.abs(parseFloat(data)) > 0.5) {
               return '<span style=\"color: red; font-weight: bold;\">' + data + '</span>';
             }
             return data;
           }"
        )
      )
    )
  ),
  rownames = FALSE,                       # Disable row names
  extensions = "Buttons"                  # Enable export buttons
) %>%
  formatStyle(
    columns = "Pdiff",
    color = styleInterval(0.99, c("black", "green")),
    fontWeight = styleEqual(0.99, "bold")  # Ensure consistent styling
  ) %>%
  formatStyle(
  columns = "deltapsi",
  color = styleInterval(c(-0.5, 0.5), c("black", "red", "black")),  # Three colors for two intervals
  fontWeight = styleInterval(c(-0.5, 0.5), c("normal", "bold", "normal"))  # Three weights for two intervals
)

write.csv(spire_alt_combination_tab,
          file = paste0("spire_alt_combination_table_", format(Sys.Date(), "%Y%m%d"), ".csv"),
          row.names = FALSE)

3.4 Gene Ontology of Filtered Splicing Events

gene_entrez_pladb <- bitr(unique(c(pladb_alt_combination_tab$GENE, pladb_exon_combination_tab$GENE, pladb_intron_combination_tab$GENE)), fromType = "ALIAS", toType = "ENTREZID", OrgDb = org.Mm.eg.db)
gene_entrez_go_pladb <- enrichGO(gene = gene_entrez_pladb$ENTREZID,
                      OrgDb = org.Mm.eg.db,
                      ont = "ALL",
                      pAdjustMethod = "BH",
                      pvalueCutoff = 0.05,
                      qvalueCutoff = 0.05)

# Polished GO dot plot
go_plot_pladb <- dotplot(
  gene_entrez_go_pladb,
  showCategory = 20
) +
  ggtitle("<b style='font-size:18px;'>GO Enrichment Analysis</b>") +
  labs(
    title = "PladB Spliced Genes",
    caption = paste0("Created by Andrés Gordo Ortiz on ", Sys.Date()),
    x = "Gene Ratio",
    y = "GO Terms"
  ) +
  theme_minimal(base_family = "sans") +
  theme(
    plot.title = element_markdown(hjust = 0.5),  # Center and style title
    axis.text.y = element_text(size = 13, color = "black"),  # Clear, professional y-axis labels
    axis.text.x = element_text(size = 12, color = "black"),  # Readable x-axis labels
    axis.title = element_text(size = 15),  # Title size enhancement
    panel.grid.major = element_line(color = "gray90"),  # Subtle gridlines
    panel.background = element_rect(fill = "white"),  # Clean background
    plot.caption = element_text(size = 10, hjust = 1, color = "gray50")  # Subtitle for the caption
  ) +
  scale_color_paletteer_c("viridis::viridis", name = "P-adjust") +  # Improve coloring for better accessibility
  scale_y_discrete(guide = guide_axis(n.dodge = 1))  # Handle overlapping labels

# FMN2DKO
gene_entrez_fmndko <- bitr(unique(c(fmndko_exon_combination_tab$GENE, fmndko_intron_combination_tab$GENE, fmndko_alt_combination_tab$GENE)), fromType = "ALIAS", toType = "ENTREZID", OrgDb = org.Mm.eg.db)
gene_entrez_go_fmndko <- enrichGO(gene = gene_entrez_fmndko$ENTREZID,
                      OrgDb = org.Mm.eg.db,
                      ont = "ALL",
                      pAdjustMethod = "BH",
                      pvalueCutoff = 0.05,
                      qvalueCutoff = 0.05)

# Polished GO dot plot
go_plot_fmndko <- dotplot(
  gene_entrez_go_fmndko,
  showCategory = 20
) +
  ggtitle("<b style='font-size:18px;'>GO Enrichment Analysis</b>") +
  labs(
    title = "FMN2DKO Spliced Genes",
    caption = paste0("Created by Andrés Gordo Ortiz on ", Sys.Date()),
    x = "Gene Ratio",
    y = "GO Terms"
  ) +
  theme_minimal(base_family = "sans") +
  theme(
    plot.title = element_markdown(hjust = 0.5),  # Center and style title
    axis.text.y = element_text(size = 13, color = "black"),  # Clear, professional y-axis labels
    axis.text.x = element_text(size = 12, color = "black"),  # Readable x-axis labels
    axis.title = element_text(size = 15),  # Title size enhancement
    panel.grid.major = element_line(color = "gray90"),  # Subtle gridlines
    panel.background = element_rect(fill = "white"),  # Clean background
    plot.caption = element_text(size = 10, hjust = 1, color = "gray50")  # Subtitle for the caption
  ) +
  scale_color_paletteer_c("viridis::viridis", name = "P-adjust") +  # Improve coloring for better accessibility
  scale_y_discrete(guide = guide_axis(n.dodge = 1))  # Handle overlapping labels

# Spire
gene_entrez_spire <- bitr(unique(c(spire_exon_combination_tab$GENE, spire_intron_combination_tab$GENE, spire_alt_combination_tab$GENE)), fromType = "ALIAS", toType = "ENTREZID", OrgDb = org.Mm.eg.db)
gene_entrez_go_spire <- enrichGO(gene = gene_entrez_spire$ENTREZID,
                      OrgDb = org.Mm.eg.db,
                      ont = "ALL",
                      pAdjustMethod = "BH",
                      pvalueCutoff = 0.05,
                      qvalueCutoff = 0.05)

# Polished GO dot plot
go_plot_spire <- dotplot(
  gene_entrez_go_spire,
  showCategory = 20
) +
  ggtitle("<b style='font-size:18px;'>GO Enrichment Analysis</b>") +
  labs(
    title = "Spire Spliced Genes",
    caption = paste0("Created by Andrés Gordo Ortiz on ", Sys.Date()),
    x = "Gene Ratio",
    y = "GO Terms"
  ) +
  theme_minimal(base_family = "sans") +
  theme(
    plot.title = element_markdown(hjust = 0.5),  # Center and style title
    axis.text.y = element_text(size = 13, color = "black"),  # Clear, professional y-axis labels
    axis.text.x = element_text(size = 12, color = "black"),  # Readable x-axis labels
    axis.title = element_text(size = 15),  # Title size enhancement
    panel.grid.major = element_line(color = "gray90"),  # Subtle gridlines
    panel.background = element_rect(fill = "white"),  # Clean background
    plot.caption = element_text(size = 10, hjust = 1, color = "gray50")  # Subtitle for the caption
  ) +
  scale_color_paletteer_c("viridis::viridis", name = "P-adjust") +  # Improve coloring for better accessibility
  scale_y_discrete(guide = guide_axis(n.dodge = 1))  # Handle overlapping labels

3.4.1 PladB

go_plot_pladb

3.4.2 FMN2DKO

go_plot_fmndko

3.4.3 Spire

go_plot_spire

3.5 Enrichment Analysis with EnrichR

3.5.1 PladB

enrichdata<-enrichr(unique(c(pladb_exon_combination_tab$GENE, pladb_intron_combination_tab$GENE, pladb_alt_combination_tab$GENE)), databases = c("GO_Biological_Process_2023","GO_Cellular_Component_2023","GO_Molecular_Function_2023"))

## Uploading data to Enrichr... Done.
##   Querying GO_Biological_Process_2023... Done.
##   Querying GO_Cellular_Component_2023... Done.
##   Querying GO_Molecular_Function_2023... Done.
## Parsing results... Done.

3.5.1.1 GO Biological Process

plotEnrich(enrichdata[[1]])

3.5.1.2 GO Cellular Component

plotEnrich(enrichdata[[2]])

3.5.1.3 GO Molecular Function

plotEnrich(enrichdata[[3]])

3.5.2 FMN2DKO

enrichdata<-enrichr(unique(c(fmndko_exon_combination_tab$GENE, fmndko_intron_combination_tab$GENE, fmndko_alt_combination_tab$GENE)), databases = c("GO_Biological_Process_2023","GO_Cellular_Component_2023","GO_Molecular_Function_2023"))

## Uploading data to Enrichr... Done.
##   Querying GO_Biological_Process_2023... Done.
##   Querying GO_Cellular_Component_2023... Done.
##   Querying GO_Molecular_Function_2023... Done.
## Parsing results... Done.

3.5.2.1 GO Biological Process

plotEnrich(enrichdata[[1]])

3.5.2.2 GO Cellular Component

plotEnrich(enrichdata[[2]])

3.5.2.3 GO Molecular Function

plotEnrich(enrichdata[[3]])

3.5.3 Spire

enrichdata<-enrichr(unique(c(spire_exon_combination_tab$GENE, spire_intron_combination_tab$GENE, spire_alt_combination_tab$GENE)), databases = c("GO_Biological_Process_2023","GO_Cellular_Component_2023","GO_Molecular_Function_2023"))

## Uploading data to Enrichr... Done.
##   Querying GO_Biological_Process_2023... Done.
##   Querying GO_Cellular_Component_2023... Done.
##   Querying GO_Molecular_Function_2023... Done.
## Parsing results... Done.

3.5.3.1 GO Biological Process

plotEnrich(enrichdata[[1]])

3.5.3.2 GO Cellular Component

plotEnrich(enrichdata[[2]])

3.5.3.3 GO Molecular Function

plotEnrich(enrichdata[[3]])

4 System Settings

sessionInfo()

## R version 4.3.1 (2023-06-16)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 22.04.4 LTS
## 
## Matrix products: default
## BLAS:   /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so.3 
## LAPACK: /usr/lib/x86_64-linux-gnu/openblas-pthread/libopenblasp-r0.3.20.so;  LAPACK version 3.10.0
## 
## locale:
##  [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
##  [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8    
##  [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8   
##  [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                 
##  [9] LC_ADDRESS=C               LC_TELEPHONE=C            
## [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       
## 
## time zone: Etc/UTC
## tzcode source: system (glibc)
## 
## attached base packages:
## [1] stats4    stats     graphics  grDevices utils     datasets  methods  
## [8] base     
## 
## other attached packages:
##  [1] clusterProfiler_4.10.1 DOSE_3.28.2            enrichplot_1.22.0     
##  [4] org.Hs.eg.db_3.18.0    org.Mm.eg.db_3.18.0    AnnotationDbi_1.64.1  
##  [7] IRanges_2.36.0         S4Vectors_0.40.2       Biobase_2.62.0        
## [10] BiocGenerics_0.48.1    patchwork_1.3.0        enrichR_3.2           
## [13] biomaRt_2.58.2         readxl_1.4.3           ggtext_0.1.2          
## [16] showtext_0.9-7         showtextdb_3.0         sysfonts_0.8.9        
## [19] ggupset_0.4.0          paletteer_1.6.0        DT_0.33               
## [22] cowplot_1.1.3          devtools_2.4.5         usethis_3.1.0         
## [25] lubridate_1.9.3        forcats_1.0.0          stringr_1.5.1         
## [28] purrr_1.0.2            readr_2.1.5            tidyr_1.3.1           
## [31] tibble_3.2.1           tidyverse_2.0.0        dplyr_1.1.4           
## [34] plotly_4.10.4          ggplot2_3.5.1          betAS_1.2.1           
## 
## loaded via a namespace (and not attached):
##   [1] splines_4.3.1           later_1.4.1             prismatic_1.1.2        
##   [4] ggplotify_0.1.2         bitops_1.0-9            filelock_1.0.3         
##   [7] cellranger_1.1.0        polyclip_1.10-7         xts_0.14.1             
##  [10] XML_3.99-0.17           lifecycle_1.0.4         lattice_0.22-6         
##  [13] MASS_7.3-60.0.1         crosstalk_1.2.1         backports_1.5.0        
##  [16] magrittr_2.0.3          sass_0.4.9              rmarkdown_2.29         
##  [19] jquerylib_0.1.4         yaml_2.3.10             remotes_2.5.0          
##  [22] rlist_0.4.6.2           httpuv_1.6.15           sessioninfo_1.2.2      
##  [25] pkgbuild_1.4.5          RColorBrewer_1.1-3      DBI_1.2.3              
##  [28] pkgload_1.4.0           zlibbioc_1.48.2         ggraph_2.2.1           
##  [31] RCurl_1.98-1.16         yulab.utils_0.1.8       WriteXLS_6.7.0         
##  [34] tweenr_2.0.3            rappdirs_0.3.3          GenomeInfoDbData_1.2.11
##  [37] ggrepel_0.9.6           tidytree_0.4.6          commonmark_1.9.2       
##  [40] codetools_0.2-20        highcharter_0.9.4       xml2_1.3.6             
##  [43] ggforce_0.4.2           tidyselect_1.2.1        aplot_0.2.3            
##  [46] farver_2.1.2            viridis_0.6.5           BiocFileCache_2.10.2   
##  [49] jsonlite_1.8.9          ellipsis_0.3.2          tidygraph_1.3.1        
##  [52] ggridges_0.5.6          tools_4.3.1             progress_1.2.3         
##  [55] treeio_1.31.0           Rcpp_1.0.13-1           glue_1.8.0             
##  [58] gridExtra_2.3           xfun_0.49               qvalue_2.34.0          
##  [61] TTR_0.24.4              GenomeInfoDb_1.38.8     withr_3.0.2            
##  [64] fastmap_1.2.0           fansi_1.0.6             digest_0.6.37          
##  [67] gridGraphics_0.5-1      timechange_0.3.0        R6_2.5.1               
##  [70] mime_0.12               colorspace_2.1-1        GO.db_3.18.0           
##  [73] markdown_1.13           RSQLite_2.3.9           utf8_1.2.4             
##  [76] generics_0.1.3          data.table_1.16.4       prettyunits_1.2.0      
##  [79] graphlayouts_1.2.1      httr_1.4.7              htmlwidgets_1.6.4      
##  [82] scatterpie_0.2.4        pkgconfig_2.0.3         gtable_0.3.6           
##  [85] blob_1.2.4              XVector_0.42.0          shadowtext_0.1.4       
##  [88] htmltools_0.5.8.1       fgsea_1.28.0            profvis_0.4.0          
##  [91] scales_1.3.0            png_0.1-8               ggfun_0.1.8            
##  [94] knitr_1.49              rstudioapi_0.17.1       tzdb_0.4.0             
##  [97] reshape2_1.4.4          rjson_0.2.23            nlme_3.1-166           
## [100] curl_6.0.1              cachem_1.1.0            zoo_1.8-12             
## [103] HDO.db_0.99.1           parallel_4.3.1          miniUI_0.1.1.1         
## [106] shinycssloaders_1.1.0   pillar_1.9.0            grid_4.3.1             
## [109] vctrs_0.6.5             urlchecker_1.0.1        promises_1.3.2         
## [112] dbplyr_2.5.0            xtable_1.8-4            evaluate_1.0.1         
## [115] cli_3.6.3               compiler_4.3.1          rlang_1.1.4            
## [118] crayon_1.5.3            labeling_0.4.3          rematch2_2.1.2         
## [121] plyr_1.8.9              fs_1.6.5                stringi_1.8.4          
## [124] viridisLite_0.4.2       BiocParallel_1.36.0     assertthat_0.2.1       
## [127] munsell_0.5.1           Biostrings_2.70.3       lazyeval_0.2.2         
## [130] colourpicker_1.3.0      Matrix_1.6-5            GOSemSim_2.28.1        
## [133] hms_1.1.3               bit64_4.5.2             KEGGREST_1.42.0        
## [136] shiny_1.9.1             gridtext_0.1.5          igraph_2.1.1           
## [139] broom_1.0.7             memoise_2.0.1           bslib_0.8.0            
## [142] thematic_0.1.6          ggtree_3.10.1           quantmod_0.4.26        
## [145] fastmatch_1.1-4         bit_4.5.0.1             gson_0.1.0             
## [148] ape_5.8